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EFFECT  OF  BLADE  STRUCTURAL  PARAMETERS  ON  HELICOPTER 


VIBRATIONAL  CHARACTERISTICS 


1.  INTRODUCTION 

This  document  represents  the  Final  Technical  Report  on  two  Army  Research 
Office  Grants,  DAAG-29-77-G-0103  and  DAAG-29-80-C-0092,  each  of  which  is  entitled 
as  per  the  title  of  this  report.  The  first  proposal  was  sent  on  8  July  1976, 
and  the  award  began  on  March  1,  1977.  This  was  to  be  a  three-year  proposal. 
However,  the  principal  investigator  spent  a  simmer  at  the  Army  Aeromechanics 
Laboratory;  and  a  no-cost  extension  was  consequently  granted  to  30  June  1980 
(a  40-month  period).  The  second  proposal  was  sent  on  1  July  1979  and  the 
grant  was  awarded  on  1  July  1980.  Again,  the  grant  period  was  initially  three 
years;  but  a  no-cost  extension  was  granted  to  31  August  1983  (a  38-month 
period).  Thus,  the  two  grants  together  span  a  period  of  78  months,  or  6  1/2 
years . 


Because  of  the  continuing  nature  of  the  second  grant,  the  Army  Research 
Office  waived  the  Final  Report  for  the  first  grant  with  the  understanding  that 
a  final  comprehensive  report  would  be  submitted  at  the  end  of  the  second  grant. 
(Reference  ARO  letter  dated  23  June  1980  DRXRO-PR,  P-14585-E,  P-1706 7-E.) 
Therefore,  this  present  document  is  the  comprehensive  report  that  covers  the 
entire  period.  This  report  treats  the  two  grants  as  a  single  unit,  and  no 
attempt  is  made  to  discriminate  work  done  before  July  1980  from  work  done 
after  July  1980. 
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2.  Statement  of  Problem 


The  major  objectives  of  our  research  project  have  been  as  follows: 

1)  Discover  the  basic  relationships  between  structural 
parameters  and  the  blade  vibrational  characteristics. 

2)  Determine  the  degree  of  sophistication  chat  is  required 
to  adequately  model  the  more  important  mechanisms. 

3)  Develop  efficient  computational  methods  suitable  to  solving 
the  rotor  trim  and  vibration  problem. 

Along  these  lines,  our  worlc  has  taken  a  "building  block"  approach  to  these 
objectives.  In  other  words,  the  statement  of  work  in  each  proposal  outlines 
small  research  tasks,  each  of  which  studies  a  particular  aspect  in  detail. 

When  two  or  more  of  these  tasks  are  completed,  the  results  are  combined 
together  to  produce  the  next  level  of  analysis.  At  each  step,  the  lessons 
learned  in  the  previous  steps  are  used  to  simplify  the  analysis  as  much  as 
possible  before  combining  with  another  analysis. 

In  the  area  of  blade  structural  vibrations,  two  separate  branches  have 
been  pursued.  One  is  the  modeling  of  blade  structural  response.  This  has 
progressed  from  rigid  blade  to  elastic  flap-lag  to  elastic  flap-lag-torsion. 

The  second  branch  has  dealt  with  problems  of  rotor-body  coupling.  Here,  we 
have  progressed  from  a  vertical  model  to  vertical-roll-pitch  to  a  complete 
9-degree-of-freedom  elastic  fuselage.  We  have  studied  both  modal  and 
finite-element  approaches. 

In  the  area  of  aeroelastic  modeling,  we  have  been  working  on  the  application 
of  simple  dynamic  stall  models  to  improve  the  estimation  of  aerodynamic  loads. 
There  are  many  sophisticated  stall  tools  that  do  not  lend  themselves  to  practical 
research  calculations,  but  we  have  limited  ourselves  to  methods  that  lend 
themselves  to  linearization. 

One  area  of  our  work  that  greatly  expanded  over  original  estimates  is 
the  study  of  rotor  trim  methods.  We  have  applied  trim  methods  over  the  entire 
spectrum  of  possible  strategies.  These  include  harmonic  balance,  numerical 
integration,  automatic  pilot,  Newton-Raphson,  periodic  shooting,  and  Floquet 
techniques. 

Finally,  because  of  the  close  relationship  between  stability,  vibration, 
and  transition  matrices,  we  have  studied  efficient  calculation  and  use  of  the 
Floquet  transition  matrix.  In  all  of  the  above  areas,  there  is  a  strong 
synergistic  relationship  among  the  areas.  Each  task  has  fed  the  other  tasks 
in  terms  of  applications  and  solutions. 
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In  addition  to  the  principal  investigator,  16  students  have  worked  on  the 
project  during  the  past  6  1/2  years.  Of  these  16  students,  12  received  advanced 
degrees  and  three  are  still  working  on  their  degrees.  Of  the  12  degree  recipients, 
four  took  positions  with  the  U.  S.  Army  or  helicopter  industry,  four  took  positions 
in  the  U.  S.  aerospace  industry,  and  four  have  continued  here  to  work  on  more 
advanced  degrees.  Below  is  a  listing  of  these  students. 


Degree 


Present  Status 


Harry  Woehrle 
Daniel  Schrage 
Abraham  Eipe 
Amir  Izadpanah 
Byung  Kim 
Timothy  Ko 
Tom  Hsu 
S-Y  Chen 
H-S  Chen 
Jon  Rogers 
Dan  Rudy 
M-S.  Huang 
Tim  Ryan 
James  O'Malley 
Swami  Karunamoorthy 
Slrajul  Iqbal 


M.S. 

12-77 

D.S. 

7-78 

D.S. 

12-79 

M.S. 

12-79 

M.S. 

5-80 

M.S. 

5-80 

D.S. 

12-80 

D.S. 

8-81 

M.S. 

12-81 

M.S. 

5-82 

M.S. 

5-83 

M.S. 

5-83 

M.S. 

8-83 

M.S. 

5-84 

U.S.  Army 
U.S.  Army 
Douglas  Aircraft 
Working  on  D.S. 
Working  on  D.S. 
Working  on  D.S. 
Sundstrand 
Kaman  Aerospace 
MTS  Corp. 

U.S.  Army 
Emerson  Space  Div. 
Working  on  D.S. 
McDonnell  Aircraft 
U.S.  Army 
Working  on  D.S. 
Unknown 
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4.  Publications  and  Reports 

Throughout  this  project,  we  have  published  the  results  of  our  research 
in  several  media.  First,  we  have  published  technical  papers  in  refereed 
journals  and  in  conference  proceedings.  Second,  we  have  produced  Master's 
and  Doctoral  theses.  Third,  we  have  published  ARO  interim  reports.  Fourth, 
we  have  participated  in  the  annual  Robert  R.  Lichten  paper  contest.  In  the 
eight  years  of  this  contest,  Washington  University  students  have  won  the 
National  Award  four  times.  The  most  recent  of  these  was  for  work  on  the 
ARO  grant.  Below  is  a  list  of  the  pertinent  publications. 
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Technical  Papers 

1.  Schrage,  D.  P.  and  Peters,  D.  A.,  "Comparison  of  the  Effect  of  Structural 
Coupling  Parameters  on  Flap-Lag  Forced  Response  and  Stability  of  a  Heli¬ 
copter  Rotor  Blade  in  Forward  Flight,"  Army  Science  Conference,  West  Point, 

New  York,  June,  1978. 

Schrage,  D.  P.  and  Peters,  D.  A.,  "Effect  of  Structural  Coupling  Parameters 
on  the  Flap-Lag  Forced  Response  of  a  Rotor  Blade  in  Forward  Flight  Using 
Floquet  Theory,"  Fourth  European  Rotorcraft  and  Powered  Lift  Aircraft 
Forum,  Stresa,  Italy,  September  1978;  Vertica,  Vol.  3,  No.  6,  June  1979. 

2.  Hsu,  T-K  and  Peters,  D.  A.,  "Coupled  Rot or /Air frame  Vibration  Analysis  by 
a  Combined  Harminic-Balance ,  Impedance-Matching  Method,"  36th  Annual 
National  Forum  of  the  American  Helicopter  Society,  May  1980,  JAHS ,  Vol.  27, 
January  1982. 

3.  Peters,  David  A.  and  Kum,  Byung  S.,  "Control  Settings  for  a  Trimmed ,  Stalled 
Rotor  by  an  Automatic  Feedback  System,”  A1AA  Dynamics  Spe  lalists'  Conference, 
Atlanta,  Georgia,  April  1981,  AIAA  Paper  No.  81-0617-CP- 

Peters,  D.  A.,  Kim,  Byung  S.,  and  Chen,  H-S,  "Calculatio.  '  Trim  Settings 
for  a  Helicopter  Rotor  by  an  Optimized  Automatic  Control.'  Journal  of 
Guidance  and  Control,  1983. 

4.  Peters,  David  A.  and  Izadpanah,  Amir,  "Helicopter  Trim  b^  .iodic  Shooting 
with  Newton-Raphson  Iteration,"  37th  Annual  Forum  of  the  American  Helicopter 
Society,  New  Orleans,  Louisiana,  May  1981,  Paper  81-23. 


Theses 


1.  Schrage,  D.  P.,  Effect  of  Structural  Parameters  on  the  Flap-Lag  Forced 
Response  of  a  Rotor  Blade  in  Forward  Flight.  Doctor  of  Science  Thesis, 
Washington  University,  May  1978. 

2.  Eipe,  Abraham,  Effect  of  Blade  Flexibility,  Structural  Parameters,  and  Trim 
Conditions  on  Rotor  Loads,  Doctor  of  Science  Thesis,  Washington  University, 
December  1979. 

3.  Izadpanah,  Amir,  Calculation  of  Helicopter  Trim  and  Air  Loads  bv  the  Method 
of  Periodic  Shooting.  Master  of  Science  Thesis,  Washington  University, 
December  1979. 

4.  Kim,  Byung,  Helicopter  Rotor  Trim  bv  an  Automatic  Feedback  System,  Master 
of  Science  Thesis,  Washington  University,  May  1980. 

5.  Ko,  Timothy,  Use  of  Tapered.  Twisted  Finite  Elements  for  Dynamic  Analysis 
of  Helicopter  Rotors.  Master  of  Science  Thesis,  Washington  University, 

May  1980. 

6.  Hsu,  Tung-Kuang,  Coupled  Rotor  Airframe  Vibration  Analysis  by  a  Combined 
Harmonic-Balance,  Impedance-Matching  Method,  Doctor  of  Science  Thesis, 
Washington  University,  August  1981. 

7.  Rogers,  Jon  P. ,  Application  of  an  Analytic  Stall  Model  to  Dynamic  Analysis 
of  Rotor  Blades,  Master  of  Science  Thesis,  Washington  University,  May  1982. 

8.  Rudy,  Daniel  J. ,  Three  Interpretations  of  a  Dynamic-Stall  Model  with 
Applications  to  Rotor  Blade  Flapping  Response,  Master  of  Science  Thesis, 
Washington  University,  May  1983. 

9.  Huang,  Ming-Sheng,  Analysis  of  Helicopter  Vibrations  with  Inplane  Degrees 
of  Freedom,  Master  of  Science  Thesis,  Washington  University,  August  1983. 


Contractor  Reports 

1.  Peters,  David  A.  and  Schrage,  Daniel  P.,  "Effect  of  Blade  Structural  Parameters 
on  the  Flap-Lag  Response  of  a  Rotor  Blade  in  Forward  Flight,"  Interim  Technical 
Report  No.  1,  ARO  Grant  DAAG-29-77-G-0103 ,  July  1978. 

2.  Peters,  D.  A.  and  Chen,  H-S ,  "Optimization  of  Auto-Pilot  Equations  for  Rapid 
Estimation  of  Helicopter  Control  Settings,  Interim  Technical  Report  No.  1, 

ARO  Grant  DAAG-29-80-C-0092 ,  November  1981. 


Invited  Lectures 


1.  Peters,  David  A.,  "Some  Mathematical  Approaches  for  Prediction  of  Rotorcraft 
Vibrations,"  presented  at  the  ARO  Rotorcraft  Vibration  Workshop,  NASA  Ames 
Research  Center,  February  1978. 

2.  Peters,  David  A.,  "Some  Aspects  of  Rotor  Trim,"  presented  at  the  ARO  Workshop 
on  Unified  Equations,  February  1979. 

3.  Peters,  David  A.,  "Helicopter  Research  at  Washington  University,  Flutter  and 
Dynamics  Council,  Denver,  Colorado,  November  1982. 

4.  Peters,  David  A.,  "The  Importance  of  Steady  and  Dynamic  Inflow  on  the 
Stability  of  Rotor-Body  Systems,"  ITR  Methodology  Assessment  Workshop, 

Ames  Research  Center,  June  1983. 


Lichten  Presentations 


1.  Hsu,  Tung-Kuang,  "Effect  of  Rotor-Body  Coupling  on  Helicopter  Vibration," 
1979. 

2.  Kim,  Byung  S.,  "Helicopter  Trim  by  an  Automatic  Pilot,"  1980. 

3.  Ko,  Timothy,  "Use  of  Tapered,  Twisted,  Finite  Elements  for  Rotor 
Blades,"  1980. 

4.  Chen,  S-Y,  "Dynamic  Analysis  of  a  Two-Bladed,  Teetering  Rotor  on  a 
Flexible  Pylon,"  1980. 

5.  Roger,  Jon,  "Application  of  an  Analytic  Stall  Model  to  Dynamic  Analysis 
of  Rotor  Blades,"  1982. 

6.  O'Malley,  James,  "Application  of  Periodic  Shooting  to  an  Existing 
Air  Loads  Program,"  1982. 

7.  Izadpanah,  Amir,  "The  Convergence  of  a  Periodic  Shooting  Algorithm  for 
Rotor  Trim,"  1982. 

8.  Rudy,  Daniel,  "Three  Interpretations  of  the  ONERA  Dynamic-Stall  Model,"  1983. 

9.  Ort,  Jack,  "Application  of  Hamilton's  Law  of  Varying  Action  to  Calculation 
of  the  Floquet  Transition  Matrix,"  1983. 

10.  Karunamoorthy ,  S.,  "Derivation  of  a  Hierarchy  of  Elastic  Blade  Equations 
for  Helicopter  Vibration  Analysis,"  1983. 

11.  Huang,  M-S,  "Analysis  of  Helicopter  Vibrations  with  Inplane  Degrees  of 
Freedom,"  1983. 


Also  presented  at 
Colorado  Springs, 


the  AIAA  Student  Conference, 
April  14-15,  1983. 


U.S. 


Air  Force  Academy, 
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5.  Summary  of  Results 


We  here  follow  the  suggestion  of  Reporting  Procedures  for  Contractors  and 
Grantees  of  the  U.S.  Army  Research  Office.  ARO  Form  18,  and  present  the  major 
results  of  our  work  in  the  form  of  reprints.  Results  not  presented  here  appear 
in  the  interim  reports  of  these  grants.  The  first  of  these  showed  that  elastic 
coupling  could  be  used  to  lower  inplane  loads,  pitch-flap  coupling  could  be  used 
to  lower  flap  loads,  and  pitch- lag  coupling  could  then  be  used  for  stability 
margins.  The  second  of  these  snowed  that  rotor  trim  could  be  obtained  by  an 
automatic  pilot  in  only  five  rotor  revolutions. 

The  reprints  given  here  represent  summaries  of  the  other  major  areas  of 
research.  The  first  paper,  by  Tom  Hsu,  summarizes  the  work  on  rotor-body 
coupling.  The  second  paper,  by  Jon  Rogers,  summarizes  the  work  on  dynamic 
stall.  Finally,  the  paper  by  O'Malley,  Izadpanah,  and  Peters  summarizes  the 
trim  studies. 


Coupled  Rotor/ Airframe  Vibration  Analysis  by  a  Combined 
Harmonic-Balance,  Impedance-Matching  Method 


T-KHsu*  D.A.  Peters 

Graduate  Student  Professor 

Washington  University,  St.  Louis,  Mo. 


A  co*  pled  rotar/sMnuoe  vibration  analysis  is  psrfonsad  by  tbs  sstebis,  of  rotas -so*  fowlofs  iMpsdssf . 
Tbs  rotor  iapstinacs  far  b  Mods*  Is  rakaJsisd  froa  tbs  psiiodk-coaffldsat  cosstloos  of  s  tiogia  Mods  lo 
forward  flight.  Thns  flapping  Modss  srs  Isdodsd,  sod  tbs  cqasdoa  b  solved  by  bsnaoak  bslsacs.  Tbs  fassiag* 
lapcdsacs.  lacMing  stractorsl  damping.  b  cskalstsd  far  3  rigid-body  sod  3  elastic  ■odes  la  pisngs,  rod.  sad 
pitch.  Tbs  rssolts  show  that  tbs  offset  of  bob  motion  os  rotor  loads  b  grsalssl  for  reistivsty  stiff  rotors,  sod  b 
sot  wsOmpproxiontsd  by  hunpsd-msii  or  portly  lasrtlsl  rotor  ■sdsb. 
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Notation 

>  cosine  and  sine  harmonics  of  F 
■slope  of  lift  curve 
■elements  of  portion  of  [2] 

■  number  of  blades 

■conventional  thrust  coefficient,  thrust /px^R4 

■  vibratory  portion  of  nondimensionai  thrust, 
roll,  and  pitch  moment  over  aa. 

■  steady  portion  of  thrust,  roil,  and  pitch  mo¬ 
ment  over  aa 

■  vector  of  harmonics  of  Cr,  CL,  CM 

=  nondimensionai  acceleration  of  gravity,  g/Q 2R 

■  plunge,  roll,  and  pitch  structural  damping 

■  fuselage  receptance 

■  hub  receptance 

■  first  flap  frequency /Q 

■  radius  of  gyration  of  fuselage  in  pitch,  roll, 
divided  by  R 

■  rotor  radius,  m 

■  fuselage  vertical  motion,  m 

■  w/R  nondimensionai  vertical  motion 
■distance  along  fuselage,  tail  to  nose,  divided  by 

R 

«  hub  vertical  displacement 
fuselage  generalized  coordinates,  Eq.  3 

-  vector  of  harmonics  of  t.  as,  ac 

-  vector  of  harmonics  of  t/,  a#,  a# 

-  rotor  impedance 

■  pitch  angle  of  hub,  fuselage,  positive  nose  up, 
rad 

■  steady  hub  pitch  angle,  rad 

■  roll  angle  of  hub,  fuselage,  positive  advancing 
side  down,  rad 
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y  -  Lock  number 

X  ■  inflow  ratio 

p  -ratio  of  distributed  beam  mass  to  total  beam 


p  -advance  ratio 

a  -rotor  solidity 

[9]  -control  derivatives 

{  0|  ■  vector  of  control  variables 

8O,0S,9C  -  collective  and  cyclic  pitch,  rad 

-  azimuth  angle,  nondimensionai  time,  ^  -  Qt 
<■>  -  swashplate  excitation  frequency  divided  by  Q 

(a>;,uy  -second  and  third  flap  frequencies,  non- 
dimensionalized  by  0 

<i„  -  frequency  of  “y”  motion  with  "x”  boundary 

condition,  divided  by  Q;  y^v.m.L— plunge, 
pitch,  roll  x = cj~ cantilevered,  free 
Q  -  rotor  speed,  rad/sec 


Introduction 

THE  concept  of  performing  a  coupled  rotor/airframe 
vibration  analysis  by  impedance  matching  goes  back  at 
least  IS  years.  Ref.  I.  That  reference  points  out  two  important 
facts.  First,  a  coupled  rotor/airframe  analysis  can  be  per¬ 
formed  in  a  rigorous  manner  by  separate  calculation  of  rotor 
and  fuselage  impedances  followed  by  a  matching  of  forces 
and  displacements  at  the  hub.  Second,  the  rotor  impedance 
need  only  be  calculated  for  a  single  blade  and  then  ap¬ 
propriately  transformed  to  apply  to  any  number  of  blades. 
Ten  years  later,  the  method  of  impedance-matching  (for¬ 
mulated  a  little  differently  from  Ref.  1)  was  used  to  illustrate 
an  important  phenomenon.1  Namely,  rotor  loads  calculated 
for  a  fixed-hub  condition  cannot  always  be  applied  as  simple 
forcing  functions  to  a  fuselage  model.  The  reason  for  this  is 
that  the  resultant  fuselage  motions  cause  the  hub  to  translate 
and  rotate  which,  in  turn,  can  alter  the  expected  loads.  This 
alteration  is  not  just  a  small  correction,  but  can  be  an  order - 
of-magnitude  change.  The  role  of  rotor  impedance  has  been 
futher  studied  by  Hohenemser*  with  very  interesting  con¬ 
clusions  that  pertain  to  fuselage  design.  In  particular,  he  notes 
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(hat  under  certain  conditions  it  may  be  desirable  to  tunc  a 
fuselage  frequency  to  the  blade  passage  frequency  in  order  to 
eliminate  hub  loads.  Hohenemser  outlines  a  method  of 
computing  the  complete  rotor  impedance  by  finite  elements 
and  transfer  matrices.  Numerical  values  of  the  effect  of  hub 
motions  on  rotor  loads  are  only  given  for  simplified  rotor 
models,  however.  Other  work  on  the  importance  of  hub 
impedance  has  appeared  in  the  literature,  for  example  Ref.  4; 
but  this  other  work  does  not  bear  directly  on  the  results 
presented  hen. 

The  interesting  conclusions  of  Refs.  2  and  3  have  raised 
some  fundamental  questions  about  coupled  rotor/body 
vibration  analyses: 

1)  If  hub  motions  have  such  a  large  effect  on  loads,  why  do 
present  methods  (which  neglect  this  effect)  show  reasonable 
correlation? 

2)  How  sophisticated  a  blade  model  is  necessary  for  a 
realistic  model  of  rotor  impedance? 

a)  Can  aerodynamic  terms  be  neglected  due  to  the  high 
frequency  ( -  4/rev)? 

b)  Can  a  constant-coefficient  or  hover  approximation  be 
used? 

c)  How  many  blade  modes  are  necessary? 

d)  Must  pitch-thrust,  roll-thrust,  and  pitch-roll  coupling 
be  included? 

3)  Which  terms  in  the  impedance  matrix  are  most 
responsible  for  the  effects  of  hub  motions  on  rotor  loads? 

It  should  be  emphasized  that  the  effects  listed  above 
(aerodynamics,  periodic  coefficients,  blade  modes,  coupling) 
have  all  been  used  previously  in  the  calculation  of  the  fixed- 
hub  loads:  but  these  terms  have  been  omitted  in  the 
calculation  of  the  rotor  impedance,  ft  is  the  purpose  of  this 
paper  to  study  the  modeling  of  rotor  impedance  and  to  answer 
the  above  questions. 

Rotor  Model 

The  rotor  model  used  here  is  that  of  Ref.  3  but  without 
reversed  flow.  The  blade  is  assumed  uniform  and  is  described 
in  terms  of  modal  coordinates  in  flapping  only.  The  periodic- 
coefficient  equations  are  solved  by  harmonic  balance  to  give  a 
solution  of  the  form 

(FI»  W(6|  +  [ZIUI  (1) 

where  (F|  is  the  vector  of  the  harmonics  of  loads,  (4|  is  the 
vector  of  controls,  and  ( z  I  is  the  vector  of  harmonics  of  hub 
motions; 


and  where  2  is  plunge,  a,  is  the  roil  angle,  and  a,  is  the  pitch 
angle  of  the  hub.  The  matrix  (4]  represents  the  rotor  loads 
with  a  fixed  hub,  and  the  matrix  |Z]  represents  the  effect  of 


hub  motions  on  rotor  toads  (i.e.,  the  impedance).  The  ef¬ 
ficiency  of  the  above  procedure  is  one  of  the  advantages  of  the 
harmonic  balance  approach.  The  calculation  of  [4]  and  1 2] 
in  Eq.  (I)  need  be  performed  for  a  single  blade  only.  Sub¬ 
sequently,  the  corresponding  matrices  for  a  b-bladed  rotor 
can  be  found  by  elimination  of  all  harmonics  (i.e.  rows  and 
columns)  that  are  not  integer  multiples  of  b.  (Complete  details 
are  in  Ref.  3). 

The  impedance  matrix  |  Z\  is  the  most  critical  item  of  the 
rotor  analysis  because  it  represents  the  crucial  effect  of  hub 
motions  on  rotor  loads.  In  Ref.  2,  this  rotor  impedance  in¬ 
cludes  the  aerodynamic  terms  for  a  rigid,  hovering  rotor  and 
the  inertial  terms  for  a  rotor  mass  lumped  at  the  hub  (i.e.  an 
infinitely  rigid  rotor),  in  the  numerical  results  of  Ref.  3, 
aerodynamic  terms  are  not  included;  but  the  inertial  terms  are 
found  with  the  assumption  that  only  one-half  of  the  rotor 
mass  is  lumped  at  the  hub  and  that  the  remaining  mass  is 
suspended  on  a  spring  of  sufficient  stiffness  to  model 
correctly  the  first  flapping  frequency.  One  of  the  purposes  of 
this  present  study  is  to  compare  the  results  calculated  under 
these  assumptions  with  results  calculated  from  the  complete, 
aeroelastic  rotor  impedance. 

It  should  be  noted  here  that  the  present  method  of 
calculation  of  rotor  impedance  has  some  experimental 
verification.  Figure  1  shows  data  for  rotor  response  due  to 
swashplate  oscillation  compared  with  harmonic-balance 
results.  Ref.  6.  The  swashplate  responses  are  sufficiently 
similar  to  shaft  responses  to  indicate  that  the  present  method 
is  reasonably  accurate  in  the  frequency  range  of  interest.  The 
peaks  in  the  curve  are  the  resonances  atp*l. 

Fuselage  Model 

The  fuselage  is  modeled  by  6  modes— a  rigid-body  and  an 
elastic  degree  of  freedom  in  each  of  the  directions  plunge, 
pitch,  and  roll.  The  plunge  model,  shown  in  Fig.  2.  is  a 
uniform  beam  with  a  lumped  mass  added  at  the  center.  The 
beam  deflection  is  expressed  in  terms  of  two  comparison 
functions  that  multiply  two  generalized  coordinates,  2  and  tf 

v(2,l)-«(t)+R(i-i4)Jf/(l)  (3) 

where  2  is  the  hub  motion,  2/  is  4  times  the  elastic  tip 
deflection,  and  t  is  the  nondimensional  distance  along  the 


Fig.  1  Compartaoa  of  theory  (with  rrnned  flow)  and  npcriiatot 
for  iwash-plale  osdllaUoa:  p  ■  2-32,  y  -  S.  p  « .71. 
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beam.  The  nondimensional  equations  of  motion  are 


‘1  p/12 

H 

0  0 

H 

'Cr 

* 

*  + 

« 

>  ■ « 

20/3  1  . 

it. 

0 

-*/- 

J »  . 

Impedance  Matching 

The  combined  rotor/ fuselage  vibrations  may  be  obtained 
by  the  matching  of  impedances  from  Eq.  (1)  with  those  of  Eq. 
(6).  This  implies  the  matching  of  harmonics  of  both  loads  { F) 
and  displacements  f  z|  at  the  interface.  Therefore,  we  have 

|Fl-{01l0)  +  (Z][ff]lF|  (7) 


where  p  is  the  ratio  of  distributed  beam  mass  to  total  beam 
mass,  <j£,  is  the  natural  frequency  of  the  fuselage  when 
cantilevered  at  the  hub,  CT  is  the  steady  thrust  required  to 
support  the  fuselage,  and  g  is  the  nondimensional  acceleration 
of  gravity  g/0 2R.  The  parameters  p  and  <j£  completely 
ifrwHiine  the  free  natural  frequency  of  the  fuselage,  J*. 

•*"5^5**-  (3) 

Structural  damping  is  included  by  multiplication  of  the 
stiffness  matrix  by  g,  to  obtain  a  damping  matrix. 

Although  Eq.  (4)  describes  a  beam  model  of  the  fuselage, 
the  equations  are  completely  analogous  to  those  for  the  much 
simpler  model  of  two  lumped  masses  connected  by  a  spring, 
as  is  used  in  Refs.  2  and  3.  In  fact,  the  fuselage  impedance  of 
Ref.  2  is  matched  by  p-  .94,  g-  .02.  Thus,  there  is  no  fun¬ 
damental  difference  between  the  vertical  fuselage  model  of 
Refs.  I  and  2  and  chat  of  this  work.  By  extension,  one  may 
also  argue  that  the  fuselage  impedance  used  here  should  be 
little  different  from  a  single-mode  impedance  obtained  from 
more  complicated,  finite-element  fuselage  models. 

The  roll  and  pitch  model  used  here  is  illustrated  in  Fig.  3  for 
the  case  of  pitch.  The  fuselage  is  considered  rigid  in  pitch  but 
connected  to  the  rotor  hub  through  a  torsional  spring.  Only 
rotational  motions  are  included,  and  hub  or  fuselage  lateral 
motion  is  reasonably  neglected.  The  resultant  equations, 
although  not  given  here,  are  exactly  analogous  to  Eq.  (4)  with 
a,  (or  ae)  taking  the  place  of  i,\  and  CL  (or  CM)  Caking  the 
place  of  Cr.  A  lumped  inertia  at  the  hub  includes  tran¬ 
smission  mass  and  hub  mass,  but  not  blade  mass  or  rotor 
moment  of  inertia.  Details  of  the  model  may  be  found  in  Ref. 
7. 

The  fuselage  impedance  for  plunge,  pitch,  and  roil  is  found 
by  substitution  of  the  appropriate  Fourier  series  for  2,  l,,  ae, 
a^,  a„  a,/,  CT,  Ct,  and  Cv  into  the  hub  equations  followed 
by  a  simple  harmonic-balance  solution. 


or 

|F)-(/-lZ][/f])-'[®ll®l  (8) 

Thus,  the  hub  loads  including  blade  motions,  (F),  are  ob¬ 
tained  from  the  fixed-hub  loads,  [0]  (01,  by  the  appropriate 
use  of  fuselage  receptance  [ATI  and  rotor  impedance  [Z\.  The 
related  hub  motions  (all  harmonics)  follow  immediately  from 
Eq.(6). 

The  efficiency  of  the  calculation  in  Eq.  (8)  should  be  em¬ 
phasized  here.  First,  we  note  that  calculation  of  rotor  im¬ 
pedance  is  by  far  the  most  difficult  portion  of  the  calculation. 
Once  [2]  and  [9]  are  calculated  for  a  given  rotor,  however, 
changes  in  fuselage  properties  [H]  can  be  made  with  only  a 
minor  effort.  Thus,  large-scale  investigations  into  the  effects 
of  fuselage  frequency  can  be  performed  with  a  very  small 
computational  effort.  When  the  coupled  rotor-fuselage 
response  is  done  in  one  computation,  on  the  other  hand,  the 
computational  effort  can  bimome  forbidding.  Another  ef¬ 
ficient  aspect  of  Eq.  (8)  is  its  potential  for  rotor  trim.  A  given 
Cr.  CL,  and  CM  (the  at  portions  of  |F))  can  be  easily 
matched  by  solving  Eq.  (8)  for  80.  8S,  and  8C  (the  first  three 
terms  of  (0|).  This  solution  involves  only  the  inversion  of  a 
3x3  matrix. 

Coupled  Response 

We  now  apply  the  foregoing  theory  to  the  calculation  of 
vibrations.  To  begin,  we  look  at  the  coupled  rotor-fuselage 
response  of  a  system  with  the  following  baseline  parameters. 

Rotor:  4  blades,  F-1.12,  «,»2.3,  u, -4.5.  y-8, 
CT,  *  .0144,  CL  »0,  CL  «0,  CM  »0,  n  m  .3,  ac  » .07. 
e».07,«w5.73,\w .032,  CT - .0058.  CT/a « .083 

Fuselage:  >>*. 37,  r^m.  14,  -  l.lSu^.  iiA-l.53, 

<ifm  -  I0.0ucm,  <i,L  -3.04, 
-4.47<id,  g,~gm  ~gL  *  .02. 


r  t*i  i 

'  (HI 

-  I*el  . 

.  Id  . 

Here,  (AfJ  and  (G1  are  receptances  (inverse  of  impedance). 

As  before,  only  integer-multiple  harmonics  of  the  blade 
number  (b.Jb.  ...  )  are  retained.  Furthermore,  higher 
harmonics  may  be  truncated  as  deemed  appropriate. 


Frequencies  with  subscript  "c”  denote  cantilevered  modes  in 
which  the  hub  degree  of  freedom  (i.e.  t,  a,,  or  ar)  is  con¬ 
strained  but  the  remainder  of  the  fuselage  is  free  to  move 
elastically.  Frequencies  with  subscript  “f”  denote  free  modes 
for  which  neither  the  hub  nor  the  fuselage  is  fixed.  For  all  of 
the  results  to  follow,  the  parameters  are  as  above  unless 
otherwise  noted. 

The  first  results  to  be  shown  are  the  4/rev  components  of 
thrust,  roll  moment,  and  pitch  moment  versus  the  un¬ 
constrained  (free)  roll  and  pitch  frequencies  and  tifm. 
Figs.  4-6.  Results  are  presented  for  p- 1.03,  1.06,  1.09,  and 
1.12.  Abo  shown  are  similar  curves,  labelled  “feedback 
neglected,”  which  give  the  fixed-hub  loads.  Several  in¬ 
teresting  characteristics  of  these  curves  should  be  noted.  First, 
the  effect  of  rotor-body  diminishes  as  ufm  (and  wa )  become 
large.  Therefore,  for  this  case,  the  rigid-body  motion  is  not 
playing  a  significant  role.  Second,  the  effect  of  coupling 
dramatically  increases  as  the  flapping  frequency  increases. 
This  indicates  that  the  rotor  inertia  plays  an  important  roll, 
but  that  the  inertia  is  isolated  from  the  4/rev  when  p  is  very 
close  to  1 .00.  Third,  the  roll  frequency  greatly  affects  the 
vertical  loads,  Cr.  This  implies  that  aerodynamic  coupling  is 
vitally  important  and  cannot  be  neglected.  Fourth,  the  roll  (or 
pitch)  moments  show  peaks  when  the  roll  frequency,  (or 
the  pitch  frequency,  <i!m )  is  near  3.7/rev;  but  they  show 
valleys  when  these  frequencies  are  at  4.0/ rev.  Thus,  as 
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Fig.  S  4/rrr  rot  lot*  a*  a  fanctloa  of  pitch  (sad  roB)  freqaeacy;  1) 
through  4)  refer  to  p- 1.03,  1.04,  1.09,  and  1.12  respectively; 


predicted  in  Ref.  3,  designing  right  on  «  blade-passage 
frequency  does  lower  hub  loads.  However,  as  we  shall  see 
later,  this  does  not  necessarily  mean  a  lowering  of  fuselage 
vibrations. 

We  now  look  at  the  effect  of  damping  on  these  results. 
Figure  7  gives  CM  versus  at  1.12  for  two  values  of 
structural  damping  g.  (Percent  critical  damping  is  roughly  j/2 
far  these  cases.)  Magnitude  and  phase  are  both  plotted. 
Several  things  are  noteworthy.  Fust,  we  see  that  with  only 
.002  damping,  the  pitch  moment  goes  to  zero,  at  up.  “4.0  as 
predicted  in  Ref.  3.  Also  noted  are  peaks  at  points  for  which 
a,*  at  iin  is  3.7.  Thus,  with  small  damping,  the  pitch-roll 
coupling  is  apparent  (two  peaks)  while  at  large  damping  it  is 
less  noticeable  (one  peak).  The  phase  behavior  is  also  in¬ 
teresting  and  shows  that  significant  phase  shifts  do  occur  near 
4/rev.  (The  phase  here  is  taken  as  the  arctangent  of  the  sin4^ 
component  divided  by  the  cos4il>  component  of  CM).  In  Fig.  8, 
we  see  the  behavior  of  the  fuselage  pitch  angle  for  the  iden¬ 
tical  range  of  frequencies.  The  curves  reveal  three  important 
behaviors.  First  of  all,  the  uncoupled  analysis  predicts  a  peak 
at  <ip,  -4.0  whereas  the  coupled  analysis  shows  the  peak  at 
d/ta-3.7.  This  is  indicative  of  the  fact  that  the  natural 
frequency  with  the  rotor  is  different  from  the  frequency 
without  the  rotor.  Second,  although  the  pitching  moment  is 
essentially  zero  at  <2^-  4.0,  gm  -  002,  the  response  at 


2.0  1.0  4.0  5. 0  0.0 

Fig.  7  4/rev  pitch  loads  as  a  fanctloa  of  pilch  (and  roll)  frequency 
■ad  damping;?- 1.12,  Jji  • 1.11  dp.. 


iifm  -4.0  is  very  high.  Therefore,  setting  CM  to  zero  by  tuning 
to  4/rev  does  not  minimize  vibrations  and  would  probably  not 
be  a  good  design  choice. 

The  physical  interpretation  of  this  is  as  follows.  First,  as  the 
natural  frequency  approaches  4.0,  the  receptance  of  the 
fuselage  goes  to  infinity  (for  zero  damping)  as  is  the  case  with 
all  oscillators.  Second,  since  H  approaches  infinity  in 
equations  (6)  and  (7),  F  must  approach  zero  (as  pointed  out  in 
Ref.  3).  However,  even  though  the  force  goes  to  zero,  the 
response  does  not  go  to  zero.  This  can  be  seen  from  Eq.  (1) 
which  shows  that  F-0  does  not  imply  2-0.  Thus  the  zero 
force  and  infinite  receptance  combine  to  create  a  finite 
response.  Also,  the  phase  angle  here  is  very  important  in  that 
fuselage  vibrations  will  be  a  combination  of  pitch  and  plunge. 
Thus,  relative  phasings  are  crucial.  A  final  look  at  the  effect 
of  damping  is  given  in  Fig.  9,  which  shows  fuselage  plunge 
amplitude  versus  pitch  frequency.  The  importance  of  dam- 
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Fig.  »  F phMfin  "  •  fmctloR  or  pilch  (and  roHI  frcqacacy 
M4  damping;  P  -  1. 12,  d/l  -  1. 1*  tifm. 


ping  and  the  strong  coupling  between  pitch  and  plunge  are 
dearly  evident.  Thus,  the  result  without  the  effect  of  hub 
motions  on  loads  is  inadequate. 

The  preceding  results  have  centered  on  the  effect  of  pitch 
and  roll  frequency.  We  now  turn  to  the  effect  of  vertical  (or 
plunge)  frequency,  J/v  Figure  10  shows  the  4/rev  Cr  asa 
function  of  d^.  The  large  differences  between  curves  with  and 
without  the  “feedback"  ofhub  motions  is  due  primarily  to  the 
effect  of  pitching  and  rail  frequency,  as  seen  previously  in 


Fig.  1#  4/rev  vortical  loads  as  a  reaction  ot  raadaga  vertical 
fraqaancy. 


trcqaaacy. 


Fig.  4.  However,  the  effect  of  elastic,  vertical  motions  is  also 
seen  to  be  significant,  as  indicated  by  the  local  decrease  in  Cr 
near  dA  -  4.0.  Two  curves  with  feedback  are  shown.  The  solid 
curve  corresponds  to  94*  of  the  fuselage  mass  distributed 
along  the  beam  and  40*  lumped  at  the  hub.  The  phase  angle 
is  fairly  constant  in  either  case.  Figure  1 1  shows  the  vertical 
hub  motion  for  the  same  parameter  range.  As  might  be  ex¬ 
pected,  the  hub  motion  is  greater  when  less  mass  is  lumped  at 
the  hub.  Although  Cr  shows  a  minimum  near  « 4.0,2 
shows  a  pronouned  peak  at  iif,  >4.0.  Figure  12  gives  the  total 
vertical  motion  along  ihe  fuselage  versus  ^  <Jt“  3  is  (he 
hub).  The  curve  for  <j,„  «4  0  is  dashed  simply  for  clarity.  The 
curves  show  that  not  only  does  the  hub  motion  peak  near 
d^  a  4.0,  but  the  tip  motion  also  peaks.  The  symmetry  of  the 
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curvet  about  “0.3  indicates,  for  this  cate,  the  relative 
dominance  of  hub  vertical  terms  over  pitch  terms,  since  both 
effects  add  to  five  the  vertical  motion  along  the  beam. 
Although  not  shown  here,  the  vibratory  pitch  and  roll 
moments  are  virtually  unchanged  by  <iy. 

So  far,  we  have  investigated  the  effects  of  the  first  blade 
flapping  frequency  and  the  effects  of  several  fuselage 
parameters.  We  now  look  at  the  effect  of  the  blade  second 
flap  frequency,  Fig.  13.  (The  third  flap  frequency,  u,,  is 
varied  as  i  .g  «,.)  The  second  flap  frequency  is  seen  to  greatly 


rtf.  14  Fnttefa  vertical  atotioa  as  a  faacttoa  of  btada  second  aatt 
Uttrtt  flap  atotte  frequency;  vj»I.S  <•,. 


affect  the  4/rev  thrust  response.  There  is  a  dramatic  rise  in 
4/rev  CT  between  u,  »  2.3  (about  the  lowest  value  physically 
possible)  and  uj  -  3.0.  The  effect  at  -  3.0  indicates  a 
response  of  the  cyclic  flapping  modes  at  (w,  +  1).  A  second 
peak  is  seen  at  u,  -4.0  and  indicates  coupling  with  a  coning 
mode.  For  >3.0,  the  effect  of  hub  motions  is  significant. 
Figure  14  gives  a  more  detailed  account  of  how  fuselage 
motions  at  various  points  along  the  beam  are  affected  by  <j,. 
The  large  dissymmetries  for  wj  >2.5  indicate  that  pitch 
moments  are  playing  a  dominant  role.  A  near  anti-node  is 
noted  at  X“.S  when  W]  is  near  3.0. 

Elements  of  Rotor  Impedance 

Now  that  we  have  taken  a  general  look  at  the  effect  of 
rotor/body  coupling  on  loads  and  vibrations,  we  are  in  a 
position  to  determine  which  terms  in  the  rotor  impedance 
matrix  are  the  primary  contributors  to  these  effects.  The  rotoi 
impedance  matrix  [Z]  is  a  rather  complicated  array,  and  we 
need  to  consider  it  in  detail  before  going  on.  In  general,  ( Z] 
represents  the  response  of  every  sine  and  cosine  harmonic  of 
CT.  CL.  and  CM  to  every  sine  and  cosine  harmonic  of  the 
motions  2,  a,,  and  a{.  For  simplicity,  let  us  momentarily 
consider  only  CT  due  to  t  for  a  b-bladed  rotor.  We  will  also 
neglect  harmonics  of  order  3b  and  higher.  This  yields 
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where  Av  are  the  elements  of  [Z] .  In  general,  dements  with 
zero  subscripts  are  not  relevant  because  higher  harmonic 
motion  has  little  effect  on  steady  loads,  and  steady  motions 
(with  the  exception  of  ar)  do  not  effect  higher  harmonic 
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The  dominant  terms  in  Eq.  (9)  are  those  terms  that  relate 
the  various  harmonics  at  the  fundamental  blade-passage 
frequency,  b/rev. 

These  terms  are  An,  At),  A„,  and  Au.  For  a  purely 
cons  taut -coefficient  system,  these  terms  would  obey  the 
equalities  A„  »A„.  A„  » -A„.  For  a  periodic-coefficient 
system,  however,  these  equalities  do  not  necessarily  hold.  The 
numerical  results  of  our  work  have  shown,  however,  that  the 
equalities  are  good  approximations  for  bz3  (rotors  with  3  or 
more  blades).  For  two-bladed  rotors,  on  the  other  hand,  the 
equalities  are  definitely  not  valid.  These  observations  could 
lead  one  to  believe  that  a  constant-coefficient  approximation 
might  be  valid  for  calculation  of  rotor  impedance  when  ba3. 
No  such  approximation  is  used  here,  however.  We  will  say 
more  about  constant-coefficients  shortly. 

In  order  to  gain  a  better  understanding  of  the  effects  of 
rotor  impedance,  pertinent  coefficients  are  defined  in  the 
following  manner.  For  example,  for  b—3  in  Eq.  (9),  the 
relationship  between  3/rev  CT  and  3/rev  Z  is  summarized  by 
the  following  coefficient. 


Cn/t,  m'Jti  (Ah  +AJ„  +Aj,  +Ajj )  (10) 

Similarly,  the  relationship  between  6/rev  Cr  and  3/rev  Z  is 
summarized  as 

Cn/Z,  •'J’AiAj,  +Ai  +AJ4I  +Al, 7  (11) 

Parallel  definitions  apply  for  other  harmonics,  for  other 
loadings  (CM  and  CL).  and  for  other  hub  motions  (a,  and  ac). 
An  interesting  aspect  of  the  coupling  between  different 
harmonics,  as  in  Eq.  (11),  is  that  such  coupling  would  be 
identically  zero  in  a  constant-coefficient  system.  The  results  to 
follow,  however,  show  that  significant  magnitudes  can  occur 
for  such  terms.  Therefore,  attempts  at  a  constant-coefficient 
approximation  for  rotor  impedance  in  forward  flight  may  not 
be  valid  in  the  presence  of  significant  higher-harmonic  loads. 
Thus,  the  tentative  conclusion  of  this  work  is  that  a  constant- 
coefficient  approximation  is  valid  for  bz3  provided  that 
2b/rev  loads  are  negligible  with  respect  to  b/rev  loads. 

Figure  13  provides  a  typical  plot  of  the  3/rev  impedance 
elements  averaged  as  in  Eq.  (10).  The  solid  lines  are  the  exact 
values  versus  advance  ratio,  and  the  dashed  lines  are  values 
calculated  from  inertial  effects  only  (aerodynamics  neglected 
as  in  a  vacuum).  The  couplings  between  vertical  and  roll-pitch 
(Cu/f„  CM ,/Zu  Cn/otjj,  Cn/a cj)  are  zero  for  the  vacuum 
calculation.  Several  interesting  observations  can  be  made 
from  this  figure  and  comparisons  with  the  results  in  Figs.  4- 
14.  First,  as  previously  mentioned,  the  vertical  frequency  does 


not  affect  roil  and  pitch  vibrations.  Therefore,  Cu/Zi  and 
Cm  Hi  must  not  significantly  contribute  to  the  role  of  rotor 
impedance  on  vibrations.  Also,  the  fact  that  roll  and  pitch 
frequencies  do  affect  thrust  indicates  that  Cn/a0  and 
Cn/ad  are  important  terms.  It  further  follows  that  the 
vacuum  approximation,  which  gives  zero  for  these  terms, 
would  not  be  adequate.  Furthermore,  this  shows  that  not  even 
the  hover  impedance,  which  also  shows  zero  coupling,  would 
be  adequate  to  model  this  important  thrust-roll-pitch 
coupling. 

There  are  two  interesting  aspects  of  Fig.  13  which  bear  on 
the  presentation  of  subsequent  results.  First,  we  note  that 
Cu/Z,  and  Cmj/Z  are  comparable  in  magnitude  to  the 
complementary  terms  Cn/au  and  Cn/acj.  Nevertheless,  we 
have  found  that  the  Z  derivatives  have  negligible  effect  on  roll 
and  pitch  whereas  the  a  derivatives  have  a  large  effect  on 
thrust.  The  reason  for  this  is  that  plunge  magnitudes  are 
generally  much  smaller  than  roll  and  pitch  magnitudes. 
Therefore,  since  the  CL  and  CM  derivatives  with  Z  have  little 
effect  or  response,  they  will  not  be  shown  in  subsequent 
graphs.  The  second  interesting  aspen  of  Fig.  13  is  the  fan  that 
Cu/ctfi m  Cmj/Oci  and  CMj/osjmCu/aej.  These  ap¬ 
proximate  equalities  have  been  found  to  hold  for  all 
parameters  and  at  all  harmonics.  Therefore,  the  two  roll- 
moment  derivatives  will  not  be  shown  in  subsequent  graphs, 
as  they  may  be  inferred  from  the  pitch-moment  derivatives. 

One  last  note  on  Fig.  IS  is  the  variation  in  curves  with  n. 
Curves  linear  with  ii  imply  an  effect  due  to  coupling  terms  of 
the  type  /tsin^  or  /icosil'.  Higher-order  curvature  represents 
the  coupling  effects  of  higher  powers  of  n  and  higher  har¬ 
monics  of  <li.  It  is  interesting  that  the  pitch  and  roll 
derivatives,  although  constant  with  n.  are  not  well  ap¬ 
proximated  by  a  vacuum  analysis.  The  error  is  207v  for 
Cmi  /orCJ  and  50V»  for  Cwj/aj. 

The  five  essential  elements  of  the  impedance  matrix, 
averaged  as  in  Eq.  (11),  are  shown  in  Fig.  16  for  a 
representative  configuration  for  the  3/rev  loads  due  to  6/ rev 
motions  and  6/ rev  loads  due  to  6/rev  motions.  The  parabolic 
characteristic  of  the  thrust/roil  couplings  indicates  the  effect 
of  higher-harmonic  periodic  coefficients.  The  total  effect  of 
these  terms  on  3/rev  response,  however,  is  dependent  on  the 
amount  of  6/rev  motion  of  Z,  cr„  and  ac.  The  counterpart 
derivatives,  6/rev  loads  due  to  3/rev  motions,  are  almost 
identical. 

Figure  17  provides  the  impedance  terms  for  4/rev  loads  due 
to  4/rev  motion.  The  curves  are  qualitatively  similar  to  those 
for  3/rev,  Fig.  13;  but  the  magnitude  of  thrust/roil  and 
thrust/pitch  coupling  is  larger  for  this  case.  Another  dif¬ 
ference  between  Figs.  13  and  17  is  that  the  CTt/Zt  derivative  is 


32 


T.  HSU 


JOURNAL  OF  THE  AMERICAN  HELICOPTER  SOCIETY 


*3.7,  tf  j  ■I.J,  a  *6.2. 


■  1  .3  .ft  .1  .1  .ft 

b)  4/ rev  WftoftwM  d—  lo  >/m  molloo. 


FI*.  M  Qm>m  of  Impedaact  mute;  p- I.1S,  «,  -3.7,  «■  -1.1, 

fotl. 


ft.  *  C  /  t,  ) 

>.  *u«*  *  <  e„  /  i 

ft.  **ut*«  *  <  t,,  /  ) 

ft.  a*Un*>  ) 

*.  —um,  ‘  1 


Flf.  II  Elements  of  Impedance  matrix,  4/rev  eomponeau  aa  a 
fanctloa  of  Lock  number,  p-  .3,  p- 1. IS  (staple  flap  mode). 


much  larger  than  it*  3/rev  counterpart  and  i*  not  well- 
approximated  by  the  vacuum  calculation.  This  increase  in 
aerodynamic  effect  is  present  despite  a  lower  Lock  number  in 
Fig.  17  than  in  Fig.  13. 

The  direct  effect  of  Lock  number  can  be  seen  in  Fig.  18. 
The  flapping  frequency  is  the  same  as  in  Fig.  17,  but  only  a 
single,  straight-line  mode  is  used.  Curves  which  approach  » 
as  7  approaches  zero  indicate  inertial  effects.  U  should  be 
pointed  out  the  Crt  refers  to  the  fourth  harmonic  of  dr/oa. 
Therefore,  y  going  to  zero  implies  "a”  going  to  zero  and 


dT/aa  going  to  infinity.  (Vacuum  curves  on  previous  plots 
assume  <?r/oo  normalized  on  a  nominal  value  of  "a”  which 
remains  constant  in  a  vacuum).  The  curves  that  remain  finite 
as  7—0  show  entirely  aerodynamic  effects. 

It  may  be  recalled  that  Figs.  4,  S,  and  6  indicate  a  strong 
effect  of  flapping  frequency  on  the  coupled  rotor/ fuselage 
response.  Figure  !9a  provides  plots  of  impedance  elements 
versus  flapping  frequency  for  a  4-bladed  rotor.  Several  in¬ 
teresting  conclusions  can  be  drawn.  First,  the  CTt/lt 
derivative  is  not  a  strong  function  of  p  in  this  frequency 
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Fit-  It  Eltmmts  of  impedance  matrix,  4/rer  components  as  a  Inaction  of  flappinf  frequency;  s'J.a*  1,0,  uj  m  2.5,  m t  ■  4.5.  a)  4/ree  com¬ 
ponents  magnitude  for  I.OOsps  1.15;  b)  4/ rev  components  phase  for  I.OOsps  1.15;  c)  4/rer  components  magnitude  for  I.DspsJ.O  (single  dap 

mode). 


range.  Therefore,  the  large  effect  of  p  on  Cr  in  Fig.  4  must  be 
attributed  to  the  coupling  terms  Cr4/a,4  and  CT4/ which 
are  dependent  on  p.  It  must  also  be  pointed  out,  however,  that 
CT4/t4  is  stongly  dependent  on  flap  frequency  at  higher 
values  of  p.  Fig.  19c.  In  fact,  for  a  single  rigid  mode  and  for  a 
uniform  blade.  CT4/t4  drops  to  zero  at  pm 2  (the  blade 
becomes  an  absorber)  and  increases  dramatically  beyond  that. 
Thus,  analyses  that  assume  a  lumped  rotor  mass  (p  -  =») 
cannot  be  expected  to  give  accurate  rotor  impedance.  Another 
interesting  aspect  of  Fig.  19a  is  the  complex  behavior  of 
several  elements  of  ( Z]  up  Os  varied.  This  indicates  complex 


dynamic  interactions  between  aerodynamic  and  inertial  terms. 
Calculations  at  higher  values  of  p  show  that  these  interactions 
are  strongly  p-dependent  even  for  p>  2.  Similar  conclusions 
hold  for  phase  angles.  Fig.  19b. 

The  final  figure  to  be  considered  here  is  the  effect  of  second 
flap  frequency  on  the  rotor  impedance.  Fig.  20.  The 
motivation  here  results  from  the  strong  effect  of  ui:  seen  in 
Fig.  13.  There  are  several  noteworthy  points  in  Fig.  20.  First, 
we  note  the  total  inadequacy  of  the  vacuum  approximation. 
The  infinite  peaks  of  that  approximation  (at  wjsJ.O, 
<i *3.0,  and  u,  *4.0)  are  simply  not  realized  in  the  presence 
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of  aerodynamic  damping.  This  is  a  companion  result  to  that 
of  Fig.  7.  The  earlier  figure  shows  that  results  without 
structural  damping  are  invalid,  and  Fig.  20  shows  that  results 
without  aerodynamic  damping  are  invalid,  even  when  w,  is 
not  directly  at  a  resonance  point.  On  the  other  hand,  despite 
the  large  effect  of  damping,  the  resonance  and  antiresonance 
effects  of  higher  modes  are  still  clearly  seen  in  the  results, 
although  diminished.  Thus,  higher  modes  play  a  large  role  in 
rotor  impedance. 

Conclusions 

Detailed  calculations  of  rotor  impedance  in  flapping  only 
have  been  calculated  for  a  wide  range  of  rotor  parameters. 
These  impedances  have  been  used  to  calculate  the  coupled 
rotor/airframe  response  of  helicopters.  The  theory,  and 
several  hundred  figures  showing  the  results,  can  be  found  in 
Ref.  7.  A  few  of  the  more  interesting  curves  have  been 


presented  here.  The  conclusions  of  these  results  are: 

1)  Hub  motions  can  greatly  affect  hub  loads  with  the  most 
dramatic  effects  being  for  relatively  stiff  rotors. 

2)  There  is  a  very  large  affect  of  pitch  and  roil  motions  on 
vertical  response,  but  little  effect  of  vertical  motion  on  pitch 
and  roil  response. 

3)  Stuctural  damping  in  the  fuselage  plays  a  dominant  role, 
with  as  little  as  1%  critical  damping  changing  the  entire 
character  of  the  fuselage  vibrations. 

4)  Higher-frequency  flap  modes  can  greatly  affect  the 
coupled  response  and  can  alter  the  relative  contributions  of 
pilch  and  plunge  motions. 

5)  Vacuum,  lumped  mass,  and  hover  approximations  of 
the  rotor  impedance  are  entirely  inadequate. 

6)  A  constant-coefficient  approximation  to  rotor  im¬ 
pedance  may  be  valid  under  certain  conditions,  but  must  be 
used  with  caution  when  higher-harmonic  loads  are  present, 
due  to  significant  cross-harmonic  coupling. 

7)  The  effects  of  Lock  number,  first  flapping  frequency, 
and  higher  flapping  frequencies  are  strongly  interdependent 
and  depend  upon  a  complicated  balance  of  inertial  and 
aerodynamic  terms. 
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Abstract 


A  dynamic  analysis  of  a  single-section  model  of  helicopter  blade  is 
performed  including  an  analytic  model  of  dynamic  stall.  This  stall  model, 
a  simplified  version  of  the  model  introduced  by  Tran  and  Petot\  character¬ 
izes  the  lift  force  on  the  blade  section  throughout  both  the  linear  and 
non-linear  regimes  of  angle  of  attack.  The  resultant,  nonlinear  blade 
equations  are  solved  by  numerical  integration  for  the  periodic,  forced 
response.  Perturbation  equations,  written  for  small  disturbances  about 
this  equilibrium,  provide  eigenvalue  and  stability  information  by  means 
of  Floquet  theory. 
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Nocation 

linear  static  lift  curve  slope,  per  degree 

linear  static  lift  curve  slope,  per  radian 

coefficients  of  the  non-linear  static  lift  curve  polynomial 

blade  element  semi-chord,  equal  to  c/2,  m 

blade  element  chord,  m 

phase  shift  parameter 

total  lift  coefficient,  lift  measured  normal  to  free  stream. 


lift  coefficient  in  linear  regime 
lift  coefficient  in  non-linear  regime 
static  lift  coefficient  in  linear  regime 
static  lift  coefficient,  approximate  expression 

actual  static  lift  coefficient 

difference  between  the  extended  linear  static  lift  coefficient 
(C  )  and  the  actual  static  lift  coefficient  (C  ) 

Z  n.  Z 

l  o 

blade  element  span,  m 

mass  moment  of  inertia  of  blade  section  about  center  of 
2  2 

rotation,  kg-ra  -  mx 

reduced  frequency  ic  =  cob/^x 

ratio  of  semi-chord  to  radial  position  =  b/x 

blade  flapping  restraint  spring,  N-m/rad,  Figure  5 
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*  mass  of  blade  section,  kg 

*  non-dimensional  blade  flapping  frequency,  per  revolution,^  i+k3/Av 

*  apparent  mass  parameter,  equation  (2) 

*  time,  sec 

■  unit  step  function 

■  forward  velocity  of  rotor  center,  m/sec 

=*  average  blade  element  velocity,  V  =»  Hx,  m/sec 
=  distance  from  center  of  rotation,  m 
3  damping  parameter,  equation  (2) 

=  blade  flapping  angle,  deg 

3 

*  Lock  number  based  on  equivalent  blade  damping,  8pabdx  /I 
=  natural  frequency  parameter,  equation  (2) 

■  parameter  relating  lift  and  airfoil  pitch  rate,  equation  (2) 

■  system  eigenvalue 

*  total  aerodynamic  angle  of  attack,  deg 

■  mean  angle  of  attack  (collective  pitch  minus  inflow  angle) ,  deg 

*  cyclic  pitch,  deg 

■  airfoil  angle  of  attack  at  which  lift  departs  from  linearity,  deg 

=  amplitude  of  airfoil  oscillation,  deg 
=  time  delay  parameter,  equation  (2) 

**  local  advance  ratio  at  blade  section,  U/ftx 

*  equivalent  advance  ratio  of  uniform  blade,  u  =  .  75jj 

■  reduced  time,  T  =  Qxt/b 

»  azimuth  angle,  i|>  ■  flt  «  kt 
-  rotor  speed,  rad/sec 

■  airfoil  frequency  of  pitch  oscillation,  rad/sec 


iv 


*  2  2 

(  ),(  )  ■  3/3t,  3  /3t  (derivatives  with  respect  to  reduced  time,  t) 

(  )  ■  periodic,  forced  response  of  (  ) 

r 

(  )  ■  perturbation  to  (  ' 


Introduction 


The  phenomenon  of  dynamic  stall  is  known  to  be  important  for  the  dynamic 
analysis  of  helicopter  rotors  in  forward  flight.  Reference  2.  Consequently,  a 
great  deal  of  effort  has  gone  into  the  development  of  analytic  methods  that 
will  predict  the  behavior  of  both  rotor  lift  and  section  pitching  moment 
during  dynamic  stall.  The  most  general  methods  available  are  based  on 
tabulated  data  in  a  three-dimensional  array  (angle  of  attack,  reduced  pitch 
rate,  reduced  pitch  acceleration)  with  an  appropriate  correction  for  local 
Mach  number,  References  3  and  4.  Other  methods  of  analysis  have  attempted  to 
replace  the  large  arrays  of  tabulated  data  with  analytic  expressions.  References 
5,  6,  and  7  outline  a  method  in  which  dynamic  stall  is  treated  as  in  impulse 
loading  that  occurs  when  stall  is  encountered  and  that  decays  following  the 
stall  event.  References  8  and  9  describe  a  method  based  on  a  dynamic  (or 
equivalent)  angle  of  attack.  The  appropriate  lift  coefficient  is  chosen  as 
the  smaller  of  two  formulas,  based  on  this  equivalent  angle.  Similar  work 
has  included  a  more  detailed  analysis  to  predict  the  onset  of  stall.  Reference  10. 
This  work  also  includes  a  time  delay  factor  to  account  for  the  finite  time 
required  to  shed  vorticity. 

Other  work  in  dynamic  stall  has  concentrated  on  the  reduction  of  the  large 
volumes  of  lift  and  moment  data  (i.e.  the  three-dimensional  tables)  to  compact 
analytic  expressions.  In  Reference  11,  dynamic  stall  data  is  represented  in 
terms  of  47  parameters  that  synthesize  the  measured  data.  In  Reference  12,  this 
same  method  is  extended  (by  the  introduction  of  an  equivalent  angle  of  attack) 
to  include  arbitrary  time  variations  in  pitching  motions.  The  equivalent  angle 
of  attack  is  based  on  Duhamel's  integral  with  the  Wagner  function  chosen  as  the 
unit  step  response  to  pitch.  The  results  indicate  good  agreement  with 


measurements . 


-2- 


In  Reference  1,  an  alternative  method  is  introduced  for  the  calculation 
of  dynamic  stall;  and  that  method  is  the  subject  of  this  paper.  The  model 
in  Reference  1  describes  the  lift  and  moment  coefficients  in  terms  of  ordinary 
differential  equations.  The  equations  are  of  first  order  in  the  linear 
regime,  and  they  are  of  third  order  in  the  nonlinear  regime.  They  include 
an  analytic  approximation  to  the  static  lift  and  pitching  moment  data  as  well 
as  12  coefficients  (6  for  lift,  6  for  moment)  that  are  chosen  based  on  experi¬ 
mental  data.  These  coefficients  are  identified  as  functions  of  angle  of  attack, 
by  parameter  identification  applied  to  +  1°  oscillations  about  each  mean  angle 
of  attack  at  various  reduced  frequencies.  The  resultant,  identified  equations 
are  used  to  simulate  lift  and  moment  for  ±6°  oscillations  about  various  mean 
angles  of  attack;  and  these  show  good  correlation  for  angles  of  attack  up  to 
23°.  Reynold's  number  effects  are  included  implicitly  (as  in  other  stall  models) 
implicitly  from  the  Reynold's  number  of  the  data  base.  Mach  number  can  be 
included  in  a  similar  manner  or  from  a  local  Mach  number  correction.  In 
Reference  13,  this  stall  model  is  applied  to  the  dynamic  response  of  a 
helicopter  blade  including  stability  and  forced  response. 

The  work  reported  in  this  paper  is  part  of  a  larger  research  effort  to 
study  the  stall  model  of  Reference  1  in  several  ways.  First,  we  wish  to 
study  the  dynamic  characteristics  of  the  stall  model  including  the  relative 
importance  of  the  parameters.  Second,  we  wish  to  study  the  implementation 
of  the  model  on  helicopter  problems.  Of  particular  interest  is  the  ease  with 
which  differential  equations  can  be  incorporated  into  conventional  rotor 
analyses.  (Such  equations  are  easily  linearized  to  obtain  conventional 
stability  information.)  Third,  we  wish  to  study  the  physical  basis  of  the 
theory  including  its  extension  to  include  unsteady  free  stream,  very  large 
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angles  of  attack,  and  plunge  (as  well  as  pitch)  motions.  This  present 
paper  is  based  on  the  work  in  Reference  14  but  addresses  only  the  lift  response 
for  a  single,  rotating  blade  section.  This  study  of  the  sectional  properties 
provides  the  foundation  for  the  continuing  work  being  done  with  the  model 
(e.g.  Reference  15)  in  which  an  entire  blade  is  considered.  Despite  the  fact 
that  the  model  in  question  has  already  been  applied  to  an  entire  rotor  blade 
in  Reference  13,  we  believe  that  much  basic  research  still  needs  to  be  done 
to  better  understand  the  model  and  to  more  consistently  apply  it  to  rotor 
problems.  We  emphasize  that  the  purpose  of  this  paper  is  not  to  assess  the 
accuracy  of  the  method  with  respect  to  other  dynamic  stall  models.  Rather, 
the  purpose  of  the  paper  is  to  study  the  basic  behavior  of  this  model,  as 
defined  in  References  1  and  13. 
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Stall  Model 

The  analytical  stall  model  presented  here  is  that  of  Reference  1. 
The  model  consists  of  three  equations  that  relate  the  lift  coefficient 
of  an  airfoil  to  its  angle  of  attack  ae  follows 


C_  +  xc,  “  kc,  +  9  +A6* 

*1  Z1  zi 


C*  +  2o^C  +  y2(l+a2)  c  -  -  Y2(l+a2)  I  AC  +  C 
z2  z2  z2  *  2 


[1C2  +  C  TT5] 


c  -  c  +  c 

z  Zl  z2 


(0 

(  3  ) 


where  C  and  C  are  the  lift  coefficients  in  the  linear  and  non-linear 
21  z2 

regions  of  angle  of  attack,  0  is  the  total  aerodynamic  angle  of  attack  of 

the  airfoil,  C  is  the  static  lift  coefficient  in  the  linear  region  of 
ZX 

angle  of  attack,  Ac  is  the  difference  between  the  extended  linear  lift 
z 

curve  (C  *  a©)  and  the  actual  static  lift  curve  (C  ),  and  C  is  the 
Zo  z  z 

*  o  _  _ 

resulting  total  lift  coefficient.  The  parameters  ,  V ,  C  are 

functions  of  blade  angle  of  attack  alone  (for  a  given  airfoil,  Reynolds 
number,  and  Mach  number)  and  must  be  determined  from  wind  tunnel  tests  by 
parameter  identification. 

The  parameters  in  equation  (1)  have  direct  physical  interpretations 
in  terms  of  classical,  unsteady  aerodynamics.  Reference  16.  The  parameter  X 
is  the  time-delay  parameter  associated  with  the  lift  deficiency  function. 

It  provides  for  changes  in  magnitude  and  phase  of  the  lift.  The  parameter  6 
represents  lift  due  to  pitch  rate.  The  parameter  -A  is  the  apparent  mass 
term.  The  parameters  in  equation  (2)  are  associated  with  the  stall 
phenomenon.  In  particular,  «  is  a  damping  parameter:  V*  is  the  frequency 
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of  Che  stall  response,  and  C  is  a  phase  shift  parameter  associated  with  the 
stall  response. 

The  above  parameters  have  been  identified  for  an  0A212  airfoil  at 

0.3  Mach  number.  Reference  17.  The  static  values  (C  ,  Ac  ,  3&C  /c>©) 

zx  z  ’  z 

are  approximated  in  Reference  1  by  polynomial  approximations.  The  other 

parameters  are  found  as  functions  of  mean  angle  attack  from  experimental 

data  at  +  1°  oscillations  about  each  mean  angle.  (The  identified  parameters 

are  later  verified  by  comparisons  with  data  at  +  6°  oscillations.)  Since 

the  tests  were  performed  over  a  range  of  reduced  frequencies,  the  resultant 

parameters  are  independent  of  reduced  frequency.  Thus,  k  *  6ob/-0.x  does  not 

enter  explicitly  into  the  equations.  It  only  enters  indirectly  (for  harmonic 

*  * 

excitation)  through  the  and  ©  terms.  Thus,  the  equations  are  applicable 
to  arbitrary  blade  motion  and  not  just  simple  harmonic  motion. 

The  above  stall  parameters  are  represented  by  simplified  functions  of  © 
that  approximate  the  true,  identified -values.  These  functions  are  given  by 
Reference  17. 

X-  0.2  (4) 

A  -  5  7T/180  (5) 

-  H"  -  ,15  11  +  l'43  <6) 

V  -  0.10  +  0.023  (©-13°)  u(0 -13°)  (7) 

cs<  -  O.  105/  p  (8) 

C  -  2  -  5.1  tan-1  {l.21(©-130)}  u(e-13°)  (9) 


The  static  lift  curve  of  the  0A212  airfoil  is  presented  in  Figure  1. 
The  curve  is  defined  in  the  linear  region  by  the  equation 


\  “  7*1  ISO  9 


9  <  0 

-  cr 


for 


(10) 
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and  in  the  non-linear  region  by  the  seventh  order  polynomidal 


C 

z 


s 


7 


1-  0 


for  6  >  9  -  10°  (11) 

-  cr 


where  the  a^  are  given  in  Reference  17  as 

a  -  1.24 
o 

al  -  0.124 

a2  -  -0.0630597  (12) 

a3  -  0.01395201 

a,  -  -0.0017390851 
4 

a5  -  0.00012451913 
a6  -  -4.6849257  x  10-6 
a?  *  7.087973  x  10-8 

As  mentioned  earlier,  the  above  stall  model (with  the  appropriate  parameters) 
has  been  applied  to  correlate  large  oscillation  wind  tunnel  test  data.  The 
results,  given  in  Reference  lfare  very  encouraging.  They  show  that  the  model 
is  indeed  capable  of  describing  the  major  dynamic  stall  characteristics  of 
airfoils. 

It  should  be  noted  here  that  the  polynomial  in  equation  (11)  diverges 

for  ©  >  30°,  as  seen  in  Figure  1.  Therefore,  for  the  work  presentee  here, 

C  is  held  constant  at  0.126  for  ©i26°. 

Equations  (1) — ( 3)  comprise  a  description  of  dynamic  stall  that  is 

amazingly  simple  to  use.  The  model  is  composed  of  differential  equations  that 

describe  lift  in  a  context  similar  to  the  manner  in  which  the  other  differential 

equations  describe  blade  motions.  C  and  C  are  alwavs  present  (no  switching 

Z1  Z2 

on  or  off  is  necessary);  but, in  the  linear  regime  ,  is  zero  so  that  has 
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no  excitation  (it  decays  to  zero  on  its  own).  Similarly,  the  C  equation 

21 

is  always  operative;  but,  in  the  unstalled  region,it  becomes  a  linear 
equation.  It  is  the  ease  of  utility  of  this  model  that  provides  the 
motivation  to  study  it. 
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Verificatioa  and  Simplification  of  Stall  Model 

The  first  step  in  the  present  application  of  the  stall  model  to  rotor 
dynamic  analysis  is  to  verify  the  validity  of  the  model  by  reproduction  of 
the  previously  published  lift  hysteresis  loops.  Thus  the  stall  model, 
equaciomO-G^has  been  transformed  into  a  computer  code  that  uses  a 
Runge-Kutta  method  of  numerical  integration  for  solution.  Figure  2  shows  the 
results  of  this  analysis  in  the  form  of  versus  B  .  The  mean  angle  of 
attack,  oscillatory  angle  of  attack,  and  blade  parameters  were  chosen  to  match 
the  examples  given  in  Reference  1  in  order  to  facilitate  the  comparison.  The 
results  of  our  computations  generally  duplicate  the  published  data  of  Reference  1 
for  mean  angles  of  attack  of  11°  or  less.  The  hysteresis  loops  tend  to  take  on 
an  ellipsoidal  appearance  when  angle-of-attack  excursions  do  not  exceed  the 
static  stall  angle  of  attack.  As  slightly  higher  mean  angles  of  attack  are 
encountered,  the  loops  tend  to  take  on  figure  eight  shapes;  and  at  mean  angles 
of  attack  well  above  stall,  the  loops  are  more  erratic.  For  mean  angles  of 
attach  between  12°  and  14°,  however,  oscillations  appear  in  the  return  portion 
of  the  present  lift  curves;  and  these  are  not  found  in  the  published  data  of 
Reference  1.  Numerical  convergence  tests  on  the  solution  have  revealed  that 
these  oscillations  do  indeed  occur  from  the  correct  numerical  solution  of 
equations  (l)-(3).  However,  if  a  relatively  course  step  size  is  used  in  the 
integration,  the  oscillations  are  effectively  filtered  out  of  the  response; 
and  the  smooth  curves  of  Reference  1  appear.  We  have  verified  this  by 
numerical  experiments  with  validated  integration  codes.  It  may  be  that  such 
an  unconscious  filtering  occurs  in  Reference  1,  thus  explaining  the  lack  of 
oscillations  in  those  results. 


The  oscillations  in  the  return  strokes  of  Figure  2  are  rather  interesting 
from  a  theoretical  standpoint  in  that  similar  oscillations  appear  in  some 
experimental  data,  such  as  that  reported  in  Reference  12.  In  other  frequency- 
response  tests,  however,  data  are  often  averaged  over  many  hysteresis  loops. 
This  averaging  (along  with  the  frequency-response  characteristics  of 
instrumentation)  sometimes  can  mask  such  oscillations.  Thus,  they  do  not 
always  appear  in  published  data.  However,  the  presence  of  these  oscillations 
in  the  data  of  Reference  12jand  in  the  model  of  Reference  1,  may  indicate 
that  the  oscillations  in  Figure  2  represent  a  true  simulation  of  a  physical 
phenomenon. 

The  next  step  in  the  research  reported  here,  is  to  study  the  stall 
model  for  possible  simplifications.  A  natural  candidate  for  simplification 

**  -k 

is  the  elimination  of  some  of  the  time  derivative  terms  involving  &  and  © 

(apparent  mass  and  angular  rate  terms).  These  terms  are  good  candidates  for 

elimination  on  at  least  two  counts.  First,  they  are  almost  always  eliminated 

in  simplified  rotor-blade  analyses;  and,  second,  their  retention  (especially  &) 

results  in  a  cumbersome  complication  in  the  state  variable  equations  for  rotor 

*  ** 

flapping.  Therefore,  it  is  useful  to  study  the  effects  of  ©  and  ©  terms  on 
the  stall  model.  In  the  first  study,  the  lift  hysteresis  loops  are  generated 
with  in  equations  (1)  and  (2)  (no  ©  or  ©  terms).  Figure  3 

shows  Cz  versus  0  for  Q  ■  6°  and  at  four  typical  mean  angles  of  attack.  The 
resulting  plots  are  not  at  all  similar  to  the  original  results  presented  in 
Figure  2.  Therefore,  the  stall  model  is  oversimplified  by  elimination  of  all 
three  parameters.  Further  investigation,  however,  shows  that  setting  only 
Che  parameter  A  to  0  produces  a  much  smaller  effect  on  the  hysteresis  loops. 

An  example  is  given  in  Figure  4  which  gives  versus  O  with  the  A  =  0, 
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k  ■  .05.  Similar  runs  have  been  made  for  0  <  <  17°  and  .05  <  k  <  .15. 

These  all  show  chat  can  be  deleted  from  the  stall  model  without  a  major 
change  in  Che  response.  One  might  argue  that  only  the  ©  should  be  removed 
from  equation  (1)  with  A.  remaining  on  the  ©  term.  As  it  turns  out, 

*  \  r 

however,  A.  is  only  12%  of  the  total  &  term  (A\>b)  so  that  either 

* 

approximation  would  be  valid,  although  including  the  ©  will  be  slightly 
more  accurate.  In  the  remainder  of  this  paper,  the  stall  model  of  Tran  and 
Petot,  equations  (l)-(3),  will  be  used  with  the  parameter  A  set  to  zero. 
This  will  be  referred  to  as  the  "simplified  stall  model." 
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Chat  the  typical  section  in  Figure  5  is  at  the  3/4  blade  radius.  Thus,  we  have 
the  definitions 


3U 

~  x 


Jj 

X  / 


where is  the  dynamic  equivalent  to  the  advance  ratio  of  a  uniform  blade, 
and  ytA  is  the  local  advance  ratio  at  the  blade  section. 

The  above  definitions  give  flapping  equations  of  the  following  form. 
(Details  are  given  in  Reference  14.) 


*  *  2  Y'  ,2  CO,  X)  /  —  i  r\2 

/3  p'/S  * -g  *  -i7-(lVl4m'',t' 


!  I  S’ 


Several  aspects  of  equation  (15)  are  noteworthy.  First,  the  equation 

results  directly  from  Newton's  laws  applied  to  the  blade  moael  in  Figure  S'  with 

o 

lift  given  by  the  lift  coefficient  multiplied  by  1/2 /o(V  +  U  sin '4')“.  Second, 
the  equation  is  written  in  terms  of  reduced  time,  "C  ,  rather  than  in  terms  of 
real  time,  t,  or  azimuth  angle,  'f'  .  The  transformation  from  X  to  either  t 
or  Y  is  easily  made,  however,  by  the  definitions  of  these  parameters. 


t 


_a  x 


— 4- 
b 


Y-  kZ-sut 


<  1  e  _(  -i 


We  have  followed  the  suggestion  of  References  (1)  and  (13)  in  using  onlv  the 
mean  portion  of  velocity , -H-X  ,  in  the  reduced  time.  (For  further  discussion  of 
this  assumption  see  Reference  15.)  Note  that  the  definition  of  k  in  equation 
(16b)  is  independent  of  any  frequency.  It  is  simply  a  measure  of  the  ratio  of 
blade  chord  to  radial  position.  A  final  comment  on  equation  (15)  is  with 
respect  to  the  parameters  "a"  and  "a".  The  parameter  a  is  the  nominal  lift- 
curve  slope  used  in  the  definition  of  '/  ,  with  one  alteration.  In  )' 


,  we  us  e 
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Rotor  Blade  Model 

The  rotor  blade  model  used  in  the  dynamics  application  of  the  simplified 
stall  model  is  presented  in  Figure  5.  The  model  consists  of  a  single  blade 
section  located  at  a  radial  position  x  from  the  center  of  rotation.  A  single 
blade  section  is  used  to  facilitate  investigation  of  the  model.  Since  any 
blade  can  be  considered  as  a  composite  of  such  sections,  this  analysis  forms 
the  foundation  for  future  work  with  an  entire  blade,  for  example  Reference  15. 

Three  virtual  forces  act  on  the  blade  element:  1)  the  centrifugal  force,  F^, 

•  • 

2)  the  D'Alembert  Force,  mx,  and  3)  the  normal  lift  force,  L,  which  is  a 
function  of  angle  of  attack,  9  ,  and  reduced  time,  1 1  .  The  blade  is  allowed 
to  flap  with  angle,  /3  ,  and  is  restrained  in  the  flapping  direction  by  a  root 
spring,  k^  .  The  angular  velocity  of  the  blade  element  about  the  hub  is  . 
Physical  dimensions  of  the  element  are  c  (the  chord),  d  (the  span),  and  b 
(the  semi-chord). 

Obviously,  the  model  depicted  in  Figure  5  is  far  from  an  actual  helicopter 
rotor  blade.  However,  please  recall  that  we  wish  to  study  the  behavior  of  the 
stall  model  for  a  typical  section.  It  follows,  that  we  would  like  the  model 
in  Figure  5  to  be  the  dynamic  equivalent  of  a  uniform  rotor  blade.  This  will 
give  the  most  realistic  response  and  the  most  valuable  comparisons.  Dynamic 
equivalence  in  hover  is  obtained  by  the  definitions 

Y  -  S/° 3  ba  x3/  Jr  ’  P2 '  J  i  +  j?rr  f  ,s"'t 

Although  this  Y*  may  look  cosmetically  different  from  the  conventional 
definition  of  Lock  number,  ft  is  actually  the  exact  dynamic  equivalent. 
Similarly  the  best  forward-flight  equivalence  is  found  under  the  assumption 
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that  the  typical  section  in  Figure  5  is  at  the  3/4  blade  radius.  Thus,  we  have 
the  definitions 


3U 

y*1  *  *4_a  x 
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3  7 


X  / 


fK  " 


where is  the  dynamic  equivalent  to  the  advance  ratio  of  a  uniform  blade, 
and  is  the  local  advance  ratio  at  the  blade  section. 

The  above  definitions  give  flapping  equations  of  the  following  form. 
(Details  are  given  in  Reference  14.) 


**  2  ^  Y  , 2  C* (e. /  _  ,  \i 


Several  aspects  of  equation  (15)  are  noteworthy.  First,  the  equation 

results  directly  from  Newton's  laws  applied  to  the  blade  model  in  Figure  S’  with 

? 

lift  given  by  the  lift  coefficient  multiplied  by  l/2/3(V  +  U  sin'f’)  .  Second, 
the  equation  is  written  in  terms  of  reduced  time, 'V  ,  rather  than  in  terms  of 
real  time,  t,  or  azimuth  angle,  vf/  .  The  transformation  from  IT  to  either  t 
or  y  is  easily  made,  however,  by  the  definitions  of  these  parameters. 

X  =  J^t  t  .k-  X  ,  Y-  Utd-cN 

We  have  followed  the  suggestion  of  References  (1)  and  (13)  in  using  only  the 
mean  portion  of  velocity ,  Ji-X  ,  in  the  reduced  time.  (For  further  discussion  of 
this  assumption  see  Reference  15.)  Note  that  the  definition  of  k  in  equation 
(16b)  is  independent  of  any  frequency.  It  is  simply  a  measure  of  the  ratio  of 
blade  chord  to  radial  position.  A  final  comment  on  equation  (15)  is  with 
respect  to  the  parameters  "a"  and  "a".  The  parameter  a  is  the  nominal  lift- 
curve  slope  used  in  the  definition  of  Y  ,  with  one  alteration.  In  ,  we  use 


"a"  as  "per  radian"  Co  conform  with  conventional  definitions.  In  the 
denominator  of  equation  (15),  however,  we  use  "a"  as  "per  degree".  This 
effectively  converts  /3  to  degrees  in  the  equation. 

The  simplified  stall  model  enters  equation  (15)  in  that 
"C  (  ©.  t  )  «  C  +  C  "is  described  by  equations  (1) — (3) ,  the  stall  model 

2  Z2 

under  investigation.  Of  importance  here  is  the  fact  that  C is  a  function  of 

*  _  * 

&  and  &  which,  in  turn,  are  functions  of  /3  and  /3  .  In  particular,  the 

total  aerodynamic  angle  of  attack  and  its  derivative  can  be  approximated  by 

the  expressions 


Q  =  ©o  +■  ©s  s  » r\  (kt)  *-  <3C  cos  C  -k  tS 
_  A/k  -t-/J 


Qi) 


9  »  0  kcoskr  -  9„ksinkt 
s  c 


-  [ 


icos(kr)  -  uB) 
1  +  psin(kT) 


(1  +  psin(kT))< 


Us  ') 


where  we  have  assumed  tan  ©  ^  ©  in  determining  the  angle  of  attack  due  to  /S  . 
(This  assumption  breaks  down,  however,  for  yM  ^  .3  as  is  discussed  in 
Reference  15.)  Equation  (17)  implies  more  than  just  a  coupling  between  the 
/3>  and  Cz  equations.  It  also  implies  a  nonlinearity  in  the  system  beyond 
the  loss  of  lift  in  stall.  To  be  precise,  we  recall  that  some  of  the 
coefficients  of  the  Cz  equation  (  &  ;  Y</  «*  ,  ^  )  depend  explicitly 

on  &  ,  equations  (6)— (9) .  Equation  (17),  therefore,  implies  more  than 
just  a  time  variation  of  these  parameters;  it  also  implies  a  nonlinear 
dependence  with  /2  . 
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One  will  also  notice  that  the  induced  flow  does  not  appear  explicitly 
in  equation  (17).  Recall,  however,  that  we  are  dealing  only  with  a  single 
section.  Thus,  any  radial  variation  of  induced  flow  does  not  enter  the 
problem.  It  follows  that  the  induced  flow  can  be  considered  as  part  of 
&0/  ©s;  and  ©c  no  loss  of  generality.  The  values  of  ©s  and  ©c 

used  in  this  present  analysis  are  based  on  approximate  trim  formulas 
developed  in  Reference  14,  which  are  based  on  a  harmonic  balance  of  the 
flapping  equations  in  the  unstalled  regime  with  and  ©o  chosen  so  as 

to  set  the  first-harmonic  flapping  to  zero.  For  small  k,  these  formulas  are 
well  approximated  by 


cm 


8  k 

3\ 


(1- 


(2.  o'} 


In  the  results  presented  here,  the  quasi-steady  version  of  equation  (20) 
is  used  (k  ■  0)  with  little  loss  of  accuracy.  It  is  emphasized,  however, 
that  equations  (19)  and  (20)  do  not  eliminate  cyclic  flapping  when  the  blade 
section  is  stalled.  A  discussion  of  trim  under  stalled  conditions  is  given 
in  Reference  15. 

Before  proceeding  to  numerical  results,  there  are  some  aspects  of  the 

theory  of  Reference  (1)  that  need  further  comment.  First,  the  model  of 

Reference  1  does  not  distinguish  between  angle  of  attack  due  to  blade  pitch 

and  angle  of  attack  due  to  vertical  velocity  components.  Both  are  treated 

identically  in  equations  (17)  and  (18).  It  is  evident  from  Reference  18, 

* 

however,  that  these  ought  to  be  treated  differently  in  the  ©  term. 
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Reference  18  convincingly  shows  that  ©  should  only  Include  the  rotation  of 

the  airfoil  with  respect  to  the  air  mass.  Thus,  only  the  geometric  pitch 

rate  should  be  included.  Furthermore,  Reference  18  shows  that  the  component 

of  along  the  blade,  /S  ,  should  also  be  included.  Thus,  a  more 

* 

consistent  &  would  be 

Q  ■*  €5>s  Jz.cosCk'fy  -  <SC  b  Sm  C-kt)  +  b/S  C<2  A 

A  second  aspect  of  the  model  of  Reference  (1)  is  that  the  theory  and  experi¬ 
ments  of  Reference  1  are  for  a  constant  free  stream,  although  the  authors  suggest 
that  the  theory  may  be  extended  to  unsteady  free  stream  by  use  of  the  average 

velocity  (in  defining  reduced  time)  and  by  implementation  of  the  unsteady  veloc- 

2 

ity  in  both  the  lift  expression,  Cz(V  +  U  sin  V  )  ,  and  in  the  angle  of 

attack.  Again,  however,  the  linear,  unsteady  aerodynamic  theory  in  Reference  19 

shows  that  this  is  not  the  best  way  to  extend  unsteady  aerodynamics  to  the  case 

of  oscillatory  free  stream.  (In  particular,  extension  to  unsteady  free  stream 

must  be  done  for  the  circulation  equations  and  not  for  lift  coefficient.) 

Nevertheless,  in  the  work  here  we  apply  the  stall  model  in  exactly  the  way 

suggested  in  Reference  1.  The  errors  introduced  by  the  above  anomalies  are 

treated  in  Reference  15  and  in  ongoing  research.  A  third  aspect  of  the  Tran 

** 

and  Petot  model  deals  with  the  apparent  mass  terms.  If  A  &  is  retained  in 
the  equations,  then  &  introduces  /3  ,  as  can  be  seen  by  inspection  of 
equation  (18).  This  necessitates  an  extra  state  variable,  ^  ,  in  the  system 
equations.  Thus,  there  is  an  advantage  to  neglect  this  term  as  we  have  done. 

As  a  final  comment  on  the  model  of  Reference  1,  we  note  that  the  model  is 
identified  for  a  nonrotating  airfoil.  One  might  expect  differences  in  the  lift 
for  rotating  and  nonrotating  airfoils  in  the  same  way  that  Loewy's  theory 
differs  from  Theodorsen's  (as  function  of  thrust  and  number  of  blades),  Reference  20. 


k. 
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However,  there  is  nothing  in  the  theory  to  imply  that  the  parameters  could  not 
be  just  as  well  identified  for  a  rotating  airfoil.  Thus,  we  proceed  with  the 
study  and  apply  the  parameters  of  equation  (4)-(9)  to  a  rotor-blade  section. 
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Tlme  History  Solutions 

The  combined  set  of  simplified-stall-model  and  blade-dynamic  equations 

form  a  fifth-order  system  (2nd  order  /3  ,  1st  order  C  ,  2nd  order  C  ). 

Z1  *  z2 

State  variables  are  introduced  into  this  system  (y^  * /3  »  y2  "  ^  ’  y3  “  Cz  » 

y.  ■  C  ,  y,  ■  C  )  and  the  resulting  equations  are  solved  by  use  of  a 
4  z2  3  z2 

predictor-corrector  msierical  integration.  The  time  history  is  begun  with 
zero  initial  conditions  and  is  continued  over  a  sufficient  number  of  cycles 
to  ensure  that  all  transients  have  decayed.  (Usually,  six  to  ten  cycles  are 
sufficient  to  obtain  convergence.)  The  forced  response  of  the  combined 
stall-model;  blade-dynamics  system  is  calculated  initially  for  the  following 
baseline  parameters 

M  -  0.30,  k  -  0.05,  V  -  6.0,  p  -  1.0,  ft-  0.20  (22a-e) 

The  choice,  k  ■  .05,  implies  b/x  ■  .05.  Later  on,  variations  will  be  made 
in  f4  and  k.  The  pitch  angle  ©o  la  specified  as  10°  so  that  the  airfoil 
will  oscillate  well  into  the  linear  and  non-linear  portions  of  the  lift  curve. 
Thus,  the  effect  of  dynamic  stall  on  the  response  of  the  airfoil  can  be  seen. 

Figure  6  presents  the  results  for  the  combined  stall-model,  blade-dynamics 
system  with  the  above  baseline  parameters.  In  Figure  6a,  the  lift  coefficient 
is  plotted  versus  ©  .  A  typical  "figure  8"  plot  is  obtained.  However,  in 
contrast  to  the  earlier  plots,  this  curve  is  complicated  by  the  fact  that  © 
includes  blade  flapping  contributions  and  is  not  purely  simple  harmonic. 

Figure  6b  presents  the  blade  flapping  response.  The  response  for  an  analysis 
without  stall  is  also  given  for  comparison.  We  notice  that  the  unstalled 
rotor  is  almost  exactly  trimmed  (no  once-per-rev  in  /S  )  whereas  the  stalled 
blade  is  not  completely  trimmed.  (There  is  a  significant  once-per-rev  component.) 
This  results  from  the  fact  that  we  are  using  approximate  trim  formulas. 
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The  average  value  of  /3  for  the  unstalled  case  Is  .13  radians,  and  the  average 
value  for  the  stalled  case  Is  .10  radians.  This  drop  in  &  can  be  attributed 
to  the  overprediction  of  lift  by  the  linearized  model.  Figure  6c  depicts  the 
variation  of  ®  with  azimuth  angle  for  both  the  stalled  and  unstalled  models. 

For  the  unstalled  case,  the  variation  in  £  is  completely  due  to  ©s  and  , 

©  ■*  v  ~  t  H  since  /3  has  little  unsteady  component. 


,  since  /3  has  little  unsteady  component. 

For  the  stalled  case,  ©*  and  remain  the  same,  but  the  cyclic  flapping 

adds  another  +  2°  to  the  angle  of  attack  yielding  ©  =  — 4>  (even  deeper 

into  stall) .  An  interesting  aspect  of  the  curves  in  Figure  6  is  the  90°  phase 

lag  between  stall  and  the  drop  in  flapping.  Under  unstalled,  trimmed  conditions, 

the  variations  in  ©  exactly  counter  the  variations  in  free  stream  to  give  a 

uniform  lift  and  uniform  /£  .  When  stall  is  introduced,  however,  there  is  a 

large  drop  in  C  at  the  maximum  ©  .  (Compare  linear  and  stall  models  in  . 
z 

Figure  6a.)  Figure  6c  shows  that  this  drop  occurs  nearly  at  vf'  ■  270°,  when 

&  is  maximum.  The  resultant  drop  in  ,  however,  occurs  about  90®  later 

at  4*  ■  360®.  This  is  the  well-known  gyroscopic  effect  that  causes  the  tip 

path  of  a  rotor  to  respond  90®  out  of  phase  with  the  lift.  Similar  analyses 

have  been  performed  over  a  range  of  advance  ratios,  yM  ,  and  nondimensional 

semi-chord,  k.  These  all  show  the  same  general  trends  as  in  Figure  6,  Reference  14. 

In  the  above  paragraphs,  we  considered  the  forced  response  of  the  coupled 

stall  model  and  rotor  model.  We  now  wish  to  turn  to  another  aspect  of  dynamic 

analysis  that  is  very  important,  the  calculation  of  transient  response 

characteristics.  First,  we  will  treat  transient  response  in  hover,  ■  0. 

Although  this  may  seem  somewhat  less  interesting  than  the  transient  response 

in  forward  flight,  it  forms  the  foundation  necessary  to  understand  the  forward 

flight  transients  (as  investigated  in  the  next  section).  To  begin,  we  give 

the  blade  an  initial  angle  of  attack  with  /S  ,  /S  ,  C  ,  C  ,  and  c  all  set 

Z1  z2  z2 

to  zero;  and  we  solve  numerically  to  investigate  the  transient  build-up  of  lift 


and  coning.  Figure  7  presents  such  results  for  three  different  step  inputs. 

First,  Figure  7a  gives  results  for  ©0  -  10®.  (Since  no  stall  occurs  for 

©  <  10®,  this  figure  is  typical  of  all  lower  values  of  &c  .)  A  typical, 

linear  transient  response  is  observed  with  a  decay  near  the  expected  value, 

e-  '  .  At  a  higher  angle  of  attack,  however,  (  ©0  *  12®  in  Figure  7b), 

the  transients  do  not  completely  decay.  Instead,  there  is  a  small  limit  cycle 

that  remains.  This  indicates  a  mild  instability  in  the  coupled  stall-flap 

equations.  When  the  blade  stalls,  begins  to  drop  which  changes  the  angle 

of  attack  causing  an  increase  in  lift.  The  flapping  angle,  /3  ,  consequently 

responds  causing  a  new  stall  cycle.  Figure  7c  shows  that  this  stall  instability 

is  more  pronounced  at  ©0  -  14®  (a  larger  limit  cycle) .  Reference  14  provides 

simultaneous  plots  of  C  and  C  which  can  be  used  with  Figure  7  to  obtain  a 

Z1  z2 

more  complete  picture  of  the  coupled  lift-rotor  instability. 
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Eigenvalue  Analysis  of  Linearized  Equations 

To  investigate  the  model  behavior  further,  both  in  hover  and  forward  fligh 
we  wish  to  write  perturbation  equations  that  will  describe  the  coupled  rotor- 
stall  model.  On  the  surface,  one  might  wonder  why  a  linearized  perturbation 
analysis  is  necessary  when  one  already  has  the  numerical  solution  of  the 
differential  equations.  Dynamic  analysis  over  the  years,  however,  has  shown 
again  and  again  that  brute-force,  numerical  solutions  of  differential  equations 
do  not  provide  the  insight  of  linearized  eigenvalue  analysis.  An  eigenvalue 
analysis  provides  specific  damping,  frequency,  and  mode  shape  data  that  is 
indispensable  for  the  understanding  of  dynamic  phenomena.  It  is  for  this 
reason  that  the  analog  and  digital  response  programs  of  the  1960's  gave  way 
to  linearized  Floquet  analysis  techniques  in  the  1970 's.  Reference  21.  It  was 
simply  impossible  to  efficiently  extract  the  necessary  insight  from  time 
histories  alone.  Thus,  dynamic  analysis  must  encompass  eigenvalue  analysis 
as  well  as  forced  response. 

For  the  model  considered  in  this  paper,  the  perturbation  equations  must 
be  written  about  a  periodic  equilibrium  that  describes  the  forced  response. 

In  other  words,  the  perturbation  equations  describe  how  the  system  behaves  when 
it  is  perturbed  away  from  the  normal  position.  (This  is  exactly  the  kind  of 
information  that  is  necessary  to  study  gust  response,  rotor-body  stability, 
or  control  derivatives  of  a  helicopter.)  In  order  to  obtain  these  equations, 
we  substitute  the  following  perturbation  expansions  into  equations  (1) ,  (2) 
and  (15)  where  subscript  F  implies  the  periodic,  forced  response  of  each 
variable  and  subscript  P  implies  the  infinitesimal  perturbation  value. 
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C.*,  -  <V*.V  •*-  (COp 

=  (C-W  +•  (C*»V 

/3  -  MV  ■*• 

©  =  (©'If  +  (&V 

As  a  result,  each  of  the  parameters  that  depend 
takes  on  expansions. 


(2Sa) 

(2lb) 

(23  c) 


(2  3cA 


upon  Q /implicitly^also 


i  -  (SW&V 
y  =  (7V*  (vO, 
<?  -  (<?V  *  (d  v 


(24^ 

(24k\ 

(2^ 


Substitution  of  equations  (23)  and  (24)  into  equations  (1),  (2),  and 
(5)  yields  a  set  of  (  )_  terms  which  by  definition  cancel  since  (  )  is  the 

r  r 

forced  solution.  The  remaining  terms  yield  a  set  of  linear  differential 
equations  in  (  )p  quantities.  Nonlinear  terms  in  (  )p  are  automatically 
eliminated  due  to  the  infinitesimal  nature  of  the  perturbation  quantities, 

(  )p.  Some  (  )p  terms  remain,  however,  as  coefficients  of  the  (  )p  quantities. 
This  introduces  additional  periodic  coefficients  in  the  equations.  A  similar 
perturbation  process  is  performed  in  Reference  13  to  obtain  linearized 
periodic-coefficient  equations.  However,  several  of  the  perturbation  terms 
have  been  inadvertently  omitted  in  Reference  13.  For  example,  in  the  C 


equation  the  perturbation  of  the  term  should  be 


rov  *Op]  ^  (evi  ~ 

h  (s'ip  £<3^p  C&)F  h(k)j>(e)p  (2£\ 
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The  (  )_,(  )-,  term  is  not  part  of  the  perturbation  equations,  the  (  )  (  )  term 

t  r  r  r 

is  negligible,  but  the  second  and  third  terms  are  the  perturbation  quantities. 
Equations  (6)  can  be  used  to  express  (&)p  in  terms  of  the  perturbation 
quantity,  (/3)p. 

(2^ 

&  =  © 

where 

.  _  —  ^  COS  jk't?)  (2 

^}P  -  Jf  /  6.VI 

•k 

The  term  (b)p(©)p  is  omitted  from  R^  and  in  equation  (11)  of  Reference  13. 
Similarly,  terms  involving  00 p  and  (C)p  are  omitted  from  the  C z  equation  in 
Reference  13.  The  correct  perturbation  equations,  given  in  Reference  14, 
contain  these  necessary  terms. 

Results  of  Eigenvalue  Analysis 

We  now  apply  the  linearized  equations  obtained  in  Reference  14  (as 

described  above)  to  study  the  behavior  of  the  coupled  stall-model,  flapping 

eigenvalues.  We  emphasize  here  that,  due  to  the  coupling  of  flapping  and 

lift,  the  eigenvalues  and  eigenvectors  of  the  system  represent  coupled  modes 

with  /3  ,  C  ,  and  C  all  participating.  Nevertheless,  one  of  these  three 
Z1  z2 

is  usually  predominant  so  that  each  mode  can  be  identified.  To  begin,  we 
examine  the  case  JU  ■  0.  For  this  case,  the  equilibrium  values  (  )p  are 
constant, which  results  in  a  set  of  constant-coefficient  equations.  The 
elgnevalue  analysis  of  these  equations  results  in  the  root  locus  plot  of 
Figure  8.  For  0°  <  <30  <  10°,  the  equations  are  completely  linear  and 
the  classical,  rigid-blade  eigenvalues  are  obtained.  The  real  part  of  the 


(ft,  ' 


c?  © 
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flapping  eigenvalue  Jj)  is  the  flap  damping^ and  the  imaginary  part  of 

the  eigenvalue  (  - -h.  \]  \  ~  ( )  is  the  damped  flap  frequency.  Eigenvalues 

also  exist  for  the  C  and  C  equations.  These  are  plotted  in  Reference  14. 

_  Z1  z2 

As  is  increased  beyond  10°,  the  frequency  and  damping  of  the  flapping  modes 

decrease  due  to  the  smaller  slope  of  C  ,  as  seen  in  Figure  1.  By  ©Q  =  12°, 

s 

the  damping  is  nearly  neutrally  stable.  It  is  this  neutral  stability  that  is 
as  a  1 t 

manifested  in  Figure  7b  A  At  *  13°,  full  stall  is  encountered  and  the 

flapping  eigenvalues  make  an  abrupt  change  of  direction.  They  quickly  become 
unstable.  It  is  this  instability  that  is  manifested  in  Figure  7c  as  a  strong 
limit  cycle.  Thus,  we  see  that  Che  time  history  behavior  is  in  agreement 
with  the  perturbation  eigenvalues. 

We  next  turn  to  the  eigenvalue  analysis  in  forward  flight.  Here,  the 

equilibrium  quantities  &S)  ,  (©)„,  (C  )_,  and  (C  )_  are  periodic.  They  are 

r  r  2^  r  Z2  F 

stored  in  files,  based  on  the  forced  response,  and  then  passed  to  the  periodic 

coefficients  in  the  Floquet  analysis.  Figure  9  presents  results  for  the 

entire  system  as  a  function  of  advance  ratio  for  an  unstalled  case.  The  real 

part  of  each  eigenvalue  is  plotted.  The  C  eigenvalue  is  the  damping  of  the 

z2 

stall  mode.  Because  the  blade  is  unstalled,  this  stall  eigenvalue  is  uncoupled 

from  the  system  and  remains  fixed  at  —  <=*  -  \0  .  There  is  a  significant 

amount  of  coupling,  however,  between  the  transient  behavior  of  C  and  that 

21 

of  /3  .  The  flapping  eigenvalue  splits  into  two  branches  at  /A  *  0.79.  This 
is  typical  ci  periodic-coefficient  systems,  e.g.  Reference  21. 

Figure  10  presents  the  effect  of  higher  collective  pitch  on  the  damping 
at  yH  *  0.25.  At  *  6°,  the  blade  begins  to  experience  significant 

stall  on  the  retreating  side.  Because  of  this,  the  flapping  and  stall  eigen¬ 
values  begin  to  interact  with  each  other  in  a  marked  way.  Both  eigenvalues  split, 
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and  Che  flap  damping  decreases  to  about  1/2  of  its  original  value  ac 
&o  m  8°.  The  type  of  damping  calculation  provided  in  Figure  10  is  possible 
because  of  Che  ease  of  implementation  of  the  stall  model  of  Reference  1  and 
because  of  its  amenability  to  linearization. 

Summary  and  Conclusions 

The  dynamic  response  of  a  single  section  of  rotor  blade  has  been  calculated 
including  an  analytic  stall  model  developed  in  Reference  1.  In  this  present 
research,  the  model  is  verified  and  simplified;  and  the  resultant  flapping- 
stall  equations  are  solved  both  by  time  history  methods  and  by  linearized, 
eigenvalue  analysis.  The  conclusions  of  these  dynamic  analyses  are: 

1)  The  simplified  stall  model  can  be  incorporated  into  a  rotor  dynamic 
analysis  by  the  addition  of  three  state  variables  at  each  section. 

(Four  are  required  for  the  complete  stall  model.)  One  such  typical 
section  is  analyzed  here. 

2)  The  coupled  blade-aerodynamic  equations  for  a  single  section  are  easily 
analyzed  for  time  history  or  for  eigenvalue  analysis,  the  latter  being 
effected  by  a  straightforward  linearization. 

3)  The  perturbation  eigenvalue  analysis  gives  results  with  direct  physical 
interpretation  with  respect  to  the  time  history  solution. 

4)  Further  research  is  necessary  to  understand  how  to  correctly  apply  the 
model  to  an  entire  blade  with  unsteady  free  stream  and  with  both  pitch 


and  plunge  motions. 
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ABSTRACT 


Most  rotor  air  loads  programs  are  trimmed  by  an  iterative 
process  with  two  steps  per  iteration.  In  the  first  step,  controls 
are  guessed  and  the  equations  are  integrated  in  time  until  all 
transients  are  decayed.  In  the  second  step,  the  controls  are 
improved  based  upon  the  difference  between  the  desired  hub  loads 
(thrust,  propulsive  force,  side  force)  and  the  computed  loads.  As 
an  alternative  to  numerical  Integration,  however,  recent  papers  nave 
suggested  a  procedure  called  periodic  shooting.  The  numerical 
shooting  procedure  can  be  used  sequentially  in  the  above,  2-step 
process;  or  it  can  be  used  in  parallel  with  the  control  strategy  as 
a  unified  trim  method. 

In  this  paper,  these  three  trim  methods  (conventional, 
sequential  shooting,  parallel  shooting)  are  applied  to  production- 
version  rotor  air  loads  programs.  The  convergence  and  efficiency  of 
the  methods  are  studied,  and  the  converged  results  are  compared  with 
wind-tunnel  data. 

1 .  Introduction 


Any  calculation  of  rotor  air  loads  requires  the  periodic 
solution  to  the  rotor  aeroelastic  equations  with  a  known  set  of 
control  settings.  Similarly,  most  dynamic  stability  calculations 
are  based  on  perturbation  equations  written  about  a  periodic 
equilibrium  position.  Therefore,  calculation  of  rotor  control 
settings  and  periodic  response  is  a  fundamental  aspect  of  rotor 
analysis. 

This  calculation  is  not  at  all  trivial,  however.  Even  when 
the  controls  are  known,  it  is  not  always  easy  to  solve  for  the 
periodic  solution.  This  is  especially  true  when  one  or  more  system 
modes  has  small  damping.  In  fact,  however,  the  rotors  controls  are 
not  known.  Instead,  what  is  known  is  the  lift  force,  propulsive 
force,  and  side  force  desired  for  a  flight  condition.  The  pilot 
controls,  therefore,  also  appear  as  unknowns  in  the  problem. 

In  general,  there  are  three  categories  of  methods  to  solve 
for  the  periodic  rotor  response.  These  are:  1)  Numerical 
Integration,  2)  Periodic  Shooting,  and  3)  Harmonic  Balance.  There 
are  also  three  categories  of  methods  for  finding  the  control 
settings.  These  are:  I)  Automatic  Pilot,  2)  Nevton-Raphson, 
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and  3)  Algebraic  Control  Equations.  Each  of  the  three  response 
methods  (1,2,3)  is  particularly  suited  for  one  of  the  three  control 
methods  (1,2,3)  in  the  sense  chat  they  are  compatible  for 
application  in  a  parallel  strategy.  For  example,  numerical 
integration  and  automatic  pilot  are  applied  in  Inference  i;  Shooting 
and  Newcon-Raphson  are  applied  in  Reference  2;  and  Harmonic  Balance 
with  Algebraic  Control  Equations  is  applied  in  Reference  3. 

Despite  this  compatibility,  however,  most  production  version 
air  loads  programs  use  numerical  integration  coupled  with  Newcon- 
Raphson  (a  rather  incompatible  combination).  The  purpose  of  this 
paper  is  to  compere  three  methods:  1)  the  conventional  numerical 
integration  with  Newcon-Raphson,  2)  the  sequential  application  of 
periodic  shooting  with  Newcon-Raphson  (without  capitalizing  on  their 
compatibility),  and  3)  the  parallel  application  of  the  two  methods. 

2.  Background 

2. 1  The  Transition  Matrix 


The  first  step  in  solution  of  a  system  of  linear  differencial 
equations  is  the  determination  of  the  transition  matrix  [$].  Given 
a  set  of  n  linear  aquations  of  the  form 

{x  }-  [A(t)]  {x}  +  b(t)  (1) 

where  A(t)  and  b(t)  are  periodic  with  period ,r  ,  the  transition 
matrix,  [$],  is  defined  such  that,  for  b(t)  ■  0, 

x(t)  =  [<D(t)]  {x(0)}  0<.t<.x  (2) 

This  further  implies  chat 


tf]  -  [AJ  [*]  (3) 

In  practice  [$]  can  be  found  by  numerical  integration  of  equation 
(1)  with  b(t)  -  0.  Tor  nonlinear  systems,  the  equation  of  state 
will  have  the  form 


{x(t)}  -  (F(x,t) } 


(4) 


It  is  often  helpful  to  linearize  these  equations  around  a  nominal  or 
periodic  equilibrium  position  {x  }.  This  solution  solves  the 
equations  P 


{ip}-  (*<*p.t)} 


(xp(0)} 


0cp(T)} 


Now,  we  write  equations  for  perturbations  about  xp(t). 

x(t)  -  x  (t)  5»(t) 

P 

where  higher  powers  of  $x  are  negligible  compared  to  5x.  Now  if 
F(x,t)  is  smooch  enough  to  have  a  Taylor  series  representation. 


(5a) 

(5b) 

(6) 

then 


(F(x.t)} 


3£i 

{F(x  ,t  )  }  ♦  [— i]  {6x(t)  } 


P  P 


Bx 


i 


(7) 
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This  leads  to  the  equations 


{<$x}-  [ — ]  tax}  (8> 

j  x-Xp 

Then  the  transition  matrix  [<p(t)J  can  be  found  from  sequential 
perturbations  of  each  element  of  (x(0)}  by  a  small  amount ,  say  ^ 
away  from  {x  (0)}.  The  resultant  perturbed  initial  conditions  can 
be  used  in  Equation  (4a) ,  and  integration  through  one  period  gives  a 
solution  (x(t)}  from  which  {6x(t)/  may  be  obtained  by 

(6x(t)}  -  {x(t)>  -  (x  (t)}  (9a) 

or  p 

(x(t)}  ■  ix  (t)}  +  (<5x(t)}  (9b) 

P 

The  transition  matrix  may  be  formed  by  dividing  the  (<5X;  columns  by  ^ 
and  assembling  them  in  [<J]  such  that 


(x(t)}  -  {x  (t)}  ♦  [4>(t)J  {6X(0)} 
P 


as  in  Equation  (7) 


2.2  Periodic  Shooting 


The  method  prescribed  here,  periodic  shooting,  utilizes  the 
transition  matrix  [  <J>]  to  find  a  periodic  solution  in  a  direct  way. 
The  first  step  in  this  procedure  (once  [<j>]  is  known)  is  to  integrate 
Equation  (1)  through  one  period  with  zero  initial  conditions  but 
with  (b(t) }  retained.  The  resultant,  non-periodic  solution  will  be 
called  {xf }. 

It  follows  from  linearity  that  the  general  solution  to 
Equation  (1)  is 


{x(t)}  -  fcf(t»  ♦  10(e)]  {x(0»  (11) 

Now  a  periodic  solution, {  x(0)  }  ■  can  be  immediately  achieved 

from  the  initial  conditions 


(  x(xl}  -{  xf(T.)  }  ♦  [0<T)  ]  {x(0) }  -  (x(0) } 

(x(0) }  (I  -  <p(x)  1  -  {xf  (x)  } 

{x(0) }  -  [I  -  o(x)]"1  ixf(t)}  (12) 

The  resultant  periodic  solution  is  obtained  from  substitution  of 
Equation  (12)  into  Equation  (11).  The  calculation  in  Equation  (12) 
is  called  "periodic  shooting"  because  the  initial  conditions  are 
"aimed"  so  as  to  hit  the  target  {x(t)}  ■  (x(0)}.  We  should  mention 
here  that  the  calculation  in  Equation  (12)  is  conceptually  identical 
(but  computationally  much  simpler)  than  the  method  described  in 
Reference  (4) . 


In  the  case  of  a  nonlinear  system.  Equation  (5a) ,  the 
procedure  is  similar  to  that  outlined  above.  For  example,  estimated 
initial  conditions,  (x^(0) },  can  be  assumed  and  an  integrated 
solution  found  (x,(t) },  that  is  not  periodic  but  is  a  first  estimate 
of  {x  }.  Thus  the  initial  conditions  can  be  modified  in  an  attempt 
to  male  (x£ }  periodic. 


58-3 


fr  (0»  S  {*^(0)}  «■  [i  (t)]'L  {x,Cr)  -  3^(0)}  (13) 

Tha  procadure  can  chan  ba  repaated  with  x  (0)  generating  a  new 
estimate  x_(t).  Thus,  the  above  algorithS  can  ba  utilized  to  find 
tha  periodic  solution  x  (t)  to  a  nonlinear  system.  It  should  be 
noted  chat  this  is  equivalent  to  a  modified  Newton-Saphson  procedure 
to  find  tha  initial  conditions  chat  will  result  in  a  periodic 
solution. 

Thus,  the  above  method  and  time-wise  integration  (until  all 
transients  decay)  stand  as  two  alternative  methods  for  the  periodic 
response.  The  third  method,  harmonic  balance,  is  not  created  in 
this  paper.  Mow,  the  complete  rotor  trim  involves  calculation  of 
control  settings  and  periodic  response.  Three  passible  means  of 
effecting  trim  are  outlined  below. 
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Figure  1.  Flow  Chart  for  Conventional  Method 


2.3  The  Conventional  Method 


This  is  e  method  which  uses  a  Newt on-Raph son  iteration 
procedure  for  convergence  on  controls  (called  control  strategy)  and 
Integrates  through  time  until  a  steady-state  solution  is  found  for 
the  given  initial  condition#.  A  flow  chart  of  this  method  is  shown 
in  Figure  1.  It  can  be  seen  that,  first,  tha  controls  are  guessed. 
Second,  the  equations  ere  integrated  In  time  until  a  periodic 
solution  is  obtalnsd  (until  all  transients  decay).  Third,  Che 
forcee  are  found.  If  they  ere  within  a  certain  error  criteria,  the 
program  stops.  If  not,  ssch  control  is  psrturbsd;  and,  for  each 
psrturbstlon,  integration  in  tlms  is  performed  until  transients 
decay.  Fourth,  a  partial-derivative  matrix  is  formed  and  new  values 
for  controls  are  found  using  a  modified  Newton-Raphson  procedure. 


y.  -i 

(8}..  ♦  ( — (  F ,  ,  -  F  .} 

old  J  desired  actual 

oqj 


(14) 
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In  some  airloads  programs  only  an  approximate  version  of 
[ljF.,/30  j~l  is  used.  Approximations  include:  1)  neglect  of 
diagonal  terms,  2)  dosed  form  approximations,  and  3)  pseudo 
inverses.  These  approximations  may  save  computation  for  eacn 
iteration,  but  at  the  possible  expense  of  requiring  more  total 
iterations. 
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Figure  2.  Flow  Chart  for  Sequential  Method 


2.4  The  Sequential  Method 

In  the  sequential  method  of  periodic  shooting,  the  right 
block  in  Figure  1  is  replaced  by  Che  shooting  algorithm  described 
previously.  This  is  depicted  in  Figure  2.  A  convergence  criteria 
must  be  applied  to  the  shooting  algorithm.  This  is  done  as  follows. 
A  solution  is  considered  to  be  converged  when  the  error  between  each 
of  the  state  variables  at  4f  ■  0  and  'F  ■  2^  is  less  than  some  chosen 
value.  Thus, every  time  the  block  diagram  calls  for  a  periodic 
solution  (i.e.,  at  every  control  perturbation),  a  new  convergence  is 
required  on  initial  conditions. 


r 


GUESS  INITIAL  CONDITIONS  ] 
SNO  CONTROLS  , 


1 

INTEGRATE  once  5 


Figure  3.  Flow  Chart  for  Parallel  Method 


2.5  Tha  Parallax  Method 


In  Chla  machod  tha  ltaracions  on  control  variables  and 
initial  conditions  in  Figure  2  are  combined  into  one  schema  chat 
iterates  simultaneously  on  controls  and  periodicity.  The  procedure 
is  similar  to  the  former  algorithm  except  chat  a  single  partial 
derivative  matrix  is  obtained  that  includes  the  changes  in  forces 
and  periodicity  with  respect  to  control  settings  and  initial 
conditions.  The  flow  chart  for  this  refined  method  is  given  in 
Figure  3.  Here,  we  have  capitalized  on  Che  fact  that  both 
strateglea  (controls  and  periodic  solution)  are  Newton-Raphson 
procedures.  Thus,  it  makes  sense  to  combine  these  into  a  single 
Newton-Raphson  scheme  with  both  controls  and  initial  conditions  as 
unknowns. 

3.  Application  to  Production  Program 
3.1  Discussion  of  Rotor  Loads  Basis 

A  test  of  the  above  described  methods  is  provided  by 
application  to  a  rotor  loads  and  performance  analysis  that  has  been 
developed  as  a  subrouting  for  use  within  the  AVRADCOM,  Applied 
Technology  Laboratory  (ATL)  V/STOL  Preliminary  Design  Program.  The 
importance  of  an  efficient  iteration  method  within  a  preliminary 
design  process  becomes  evident  as  tha  analysis  is  permitted  to  allow 
more  and  more  design  variables  to  be  considered.  The  basis  for  Che 
applied  rotor  analysis  is  documented  in  the  cited  References  5,  6. 
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the  basic  equations  are  for  a  rigid  hinged  blade  with  hinge  offset. 
Only  flapping  dynamics  are  considered.  Calculation  of  the  rotor 
loads  requires  the  airfoil  section  lift  and  drag  characteristics  as 
well  as  the  resultant  velocity.  The  airfoil  section  characteristics 
are  provided  for  section  angles  of  attack  from  -180°  to  180°  for 
Mach  numbers  from  0.  to  1.0.  The  use  of  basic  steady  airfoil  lift 
and  drag  measured  as  a  result  of  2-  dimensional  transonic  wind 

tunnel  tests  is  done  with  confidence  along  most  of  the  rotor  blade 
span.  However,  three  separate  adjustments  are  required  to  account 
for  air  flow  and  blade  motion  which  can  become  significant  depending 
upon  the  rotor  operating  regime.  The  first  of  these  adjustments  is 
a  so-called  tip  relief  model,  derived  in  Reference  7,  which  accounts 
for  the  reduced  compressibility  existing  in  the  3-dimensional  flow 
near  the  tip  of  a  lifting  surface.  The  tip  relief  model  is  based 
upon  the  potential  representation  of  the  thickness  effect  of  an 
airfoil  by  a  source-sink  distribution.  The  thickness  effect  can  be 
thought  of  as  a  qualitative  explanation  for  tip  relief  in  the  3ense 
that  2-dlmenslonal  flow  requires  greater  displacement  in  a 
perpendicular  direction  than  3-dimensional  flow  about  a  finite  tip. 
Therefore,  there  results  a  relief  in  the  flow  about  the  tip  as 
compared  to  the  2 -dimensional  flow.  The  potential  function  is 
formulated  for  a  finite  wing  by  subtracting  the  functions  for 
complementary  wings  on  both  sides  from  the  function  for  an  infinite 
wing  (2-dlmenslonal  airfoil) .  Formulation  in  this  manner  relates 
the  velocity  on  the  2-dimsnsional  airfoil  to  that  on  the  finite 
wing. 


The  second  adjustment  to  the  2-dimensional  airfoil  data  is 
intended  to  account  for  the  radial  flow  conditions  that  exist  on  a 
rotor  blade  due  to  its  yawed  position  present  for  much  of  the 
azimuthal  circuit.  The  significant  features  of  this  method 
(Reference  8)  include  an  estimate  for  the  increased  skin  friction 
drag  due  to  the  use  of  the  resultant  velocity  acting  at  a  yaw  angle 
to  the  blade  element  and  a  stall  delay  due  to  an  Increased  lift 
capability  evidenced  in  yawed  flow  experiments  on  various  wings. 

The  third  adjustment  to  the  basic  wind  tunnel  tested  airfoil 
data  la  an  approximation  of  the  stall  hysteresis  with  lift  overshoot 
that  occurs  as  a  result  of  an  airfoil  oscillating  near  stall.  The 
cyclic  pitch  variation  required  by  a  conventional  rotor  system 
causes  this  unsteady  airfoil  characteristic  to  have  a  significant 
effect  upon  calculations  when  the  operating  condition  allows 
appreciable  stall.  The  formulation,  detailed  in  Reference  8,  is 
based  upon  tests  of  four  airfoil  sections  from  6Z  to  12%  thick. 
Derived  from  these  teats  are  the  stall  delay  angles  as  a  function  of 
a  dimensionless  parameter,  /1C&/2V  j  (analogous  to  the  reduced 
frequency  parameter)  where  C  -  blade  chord  and  V  -  local  velocity. 
Linear  functions  have  been  developed  for  a  stall  delay  parameter 
which  depend  on  the  airfoil  thickness,  Mach  number,  2V| 
parameter,  and  whether  it  is  lift  or  moment  stall  which  is  being 
examined.  The  moment  stall  formulation  is  used  to  determine  the 
unsteady  drag  coefficient.  Reference  9  shows  this  to  be  a  good 
approximation. 

The  non-dimensional  integral  expressions  for  the  three  rotor 
forces  (thrust  along  rotor  shaft,  and  propulsive  force  and  side 
force,  perpendicular  to  each  ocher  and  the  rotor  shaft)  are  derived 
from  the  resolution  of  the  airfoil  force  coefficients  as  they  vary 
along  the  rotor  blade.  The  integral  spans  the  distance  from  the 
root  cutout  r  to  the  tip  (1).  Tip  losses,  or  the  approximation  of 
Che  loss  of  lift  due  to  the  finite  blade,  are  approximated  by 
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setting  lift  -  0.0  at_r  ■  1  ana  assuming  a  ilnssr  variation  in  che 
lift  from  f  »  .97  to  c  »  i.  Ths  drag  fores  cosfficisne  used  ac 
r  *  l  is  chac  which  has  bssn  "aieuiaesd  as  a  result  of  applying  che 
abovs  summarized  dp  rallsr  method  ae  lift  ■  0.0. 

3.2  Application  of  learacion  Machoas 

Tha  application  of  tha  procadura  suosarized  above  requires  an 
iteration  method  to  solve  for  cha  required  rotor  forces  and  cha 
accoir  anying  staady-scace  rotor  blade  motion.  Specifically,  cha  ATL 
V/STuu  Preliminary  Design  Program  requires  rotor  tOTqua  and  cip 
path-plans  inclination  when  given  Che  rotor  forces.  Tha  iteration 
method  must  provide  convergence  of  magnitude  and  direction  upon  che 
resultant  of  rotor  lift*  propulsive  force  and  aide  force  (Figure  4) . 
This  is  done  by  adjusting  collective  pitch  ( Q  ) ,  longitudinal  cyclic 
angle  (6  ),  and  lateral  eyede  (6  ).  Steady  state  blade  flapping 
magnltudi  and  velocity  oust  be  attained.  Three  methods  of  Iteration 
have  been  applied  to  investigate  che  relative  efficiency  of  each  in 
achieving  convergence  through  Cha  variation  of  che  five  variables. 


V/OMH  MAT  ALPHAS  CLR/3  CXR/3 

0.300  0.778  -1000  0.0718  0.003380 


-3  1  ■  r— — — —  i-  »'  1  r "  -  i  —I  —  .  .  —  I  ■  i 

a  a  a  a  a  ■■  io  it  la 

ROTOR  REVOLUTIONS 


Figure  5.  Step  by  Seep  Conventional  Method 


1)  Step  by  step  conventional  method:  This  first  method 
(Figure  3)  steps  through  each  one  of  the  variables,  insuring 
convergence  within  a  specified  tolerance  before  proceeding  to  the 
next  variable.  Steady  State  flapping  is  calculated  first  and  then 
©  is  incremented  toward  the  resultant  vector  magnitude  convergence. 
Total  force  is  not  integrated  until  a  steady  state  flapping  is 
achieved.  When  the  force  magnitude  is  converged, 0  is  incremented 
toward  the  vector  direction  required.  The  vector  mlgnltude  is  then 
checked  and  then  reiterated  until  it  again  is  converged.  This 
procedure  is  repeated  until  both  magnitude  and  direction  are 
correct.  At  this  point, &  is  incremented  and  when  the  results  of 
this  perturbation  are  available,  tests  are  made  to  check  the 
previously  achieved  convergence  on©  and  &  .  If  this  test  shows 
non-convergence  the  procedure  is  repeated  from  the  point  of 
non-convergence.  Extrapolation  and  interpolation  are  accomplished 
in  small  enough  linear  steps  so  as  to  approximate  the 
non-linearities  of  the  problem.  When  enough  consistent 
perturbations  have  been  accomplished,  the  step  by  step  procedure  is 
deviated  upon,  in  that  when  one  control  is  incremented,  enough  is 
known  about  the  sensitivities  so  that  the  other  controls  can  be 
changed  at  the  same  time.  Thus  the  off-diagonal  terms  (coupling 
terms)  of  the  Inverse  matrix  can  be  included.  Upon  changing  €*  : 
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where  1  and  2  are  from  magnitude  or  direction  converged  conditions. 


v/oua  JOfT  ALPHAS  CLR/3  CXR/3 

0300  3.770  -3.000  0.0719  3.000300 


Figure  6.  Newton-Raphson  Sequential  Method 


2)  Sequential  Newton-Raphson:  This  method  (Figure  6) 
requires  the  establishment  o£  a  matrix  of  slopes  representing  the 
sensitivity  of  the  rotor  forces  to  the  isolated  perturbation  of  each 
of  the  control  variables.  This  matrix  Is  applied,  through 
Inversion,  to  the  Newton-Raphson  equation  to  achieve  simultaneous 
convergence.  The  method  is  termed  Sequential  because  each 
perturbation  requires  first  the  convergence  of  flapping  displacement 
and  velocity.  Flapping  convergence  Is  achieved  through  periodic 
shooting  with  the  first  perturbation  being  the  value  found  from 
numerical  Integration.  If  an  accurate  first  estimate  is  possible 
for  Che  flap  motion,  numerical  integration  will  yield  a  periodic 
solution  within  two  revolutions  for  a  practical  articulated  rotor. 

For  the  two  variables  included  in  this  problem,  sequential  periodic 
shooting  would  require  four  rotor  revolutions  Co  establish  the 
required  matrix  and  then  another  (minimum)  to  converge,  for  a  total 
of  five  revolutions.  Because  convergence  is  tested  for  every 
perturbation,  the  Newton-Raphson  sequential  integration  (conventional 
method)  is  superior,  in  this  application,  to  Che  Newton-Raphson 
sequential  shooting  technique.  Sequential  shooting  would  be 
advantageous  for  more  blade  motion  degrees  of  freedom.  When  all  the 
control  variables  are  perturbed  and  the  simultaneous  solution  of  all 
Che  variables  does  not  result  in  convergence,  two  variations  of 
matrix  update  are  possible.  The  first  variation  checks  to  see  which 
parameter  is  furthest  from  convergence  and  then  allows  a 
perturbation  of  this  single  variable  in  order  to  update  only  the  one 
affected  matrix  column.  A  new  es time  -on  is  then  made  for  all  the 
variables  and  convergence  la  retested.  The  second  variation 
requires  perturbation  of  each  of  the  variables  when  convergence  is 
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not  achieved,  thus  resulting  in  a  completely  updated  matrix.  This 
second  variation  is  the  procedure  for  which  results  are  presented. 


V/OMK  MAT  ALPHAS  CLR/S  CXR/S 

0.300  0.779  -3.000  0.0713  0.003300 


Figure  7.  Newton-Raphson  Parallel  Shooting 


3)  Parallel  Newton-Raphson  shooting:  This  method  (Figure  7) 
extends  Che  sequential  method  to  include  in  the  Newton-Raphson 
formulation,  Che  periodic  blade  motion  variables.  Isolated 
perturbation  of  each  of  the  variables  is  used  to  construct  a 
combined  matrix  of  slopes  representing  the  sensitivity  of  both  the 
forces  and  periodic  blade  motion.  The  term  "parallel"  then,  refers 
to  Che  fact  that  the  periodic  blade  motion  is  being  iterated  upon  at 
the  same  time  as  the  integrated  force  magnitude  and  direction. 

3.3  Method  of  Application 

Each  of  the  above  three  methods  has  been  used,  on  an  equal 
basis,  in  conjunction  with  the  rotor  loads  and  performance  method 
summarized  above.  An  equal  basis  of  comparison  is  assured  by  the 
use  of  Che  same  estimates  for  the  starting  values  of  the  control 
variables.  For  each  series  of  calculations  (each  rotor  shaft  angle 
at  a  particular  advance  ratio)  the  first  point  uses  the  estimates 
for  control  angles  and  blade  motion  based  on  a  closed  form  solution. 
The  subsequent  points  use  this  same  closed  form  solution  for  the 
angles,  but  Che  estimate  is  modified  based  on  the  differences 
between  the  initial  values  and  the  converged  values  for  the  previous 
point.  Improvements  to  this  scheme  would  include  extrapolating  the 
converged  controls  (angles  and  blade  motion)  based  on  the  previous 
two  values.  The  convergence  criteria  is  the  same  for  all  cases:  1% 
of  resultant  force  magnitude;  1%  of  resultant  force  direction;  10% 
of  rotor  side  force;  1%  of  blade  flapping  angle  and  velocity  (except 
for  small  angle;  the  tolerances  for  which  is  .001  radians).  These 
convergence  criteria  are  small  enough  such  that  a  consistent  set  of 
data  can  be  calculated.  Overall  rotor  performance  is  relatively 
insensitive  to  side  force,  so  the  larger  tolerance  is  acceptable. 
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Each  perturbation  includes  a  tolerance  test  on  every  required 
variable,  so  eacn  time  a  perturbation  is  required,  the  control  is 
perturbed  in  the  direction  toward  convergence.  Each  time  the 
controls  are  calculated  as  a  result  of  the  complete  matrix  update, 
they  are  saved  to  be  used  with  their  previous  counterparts  to 
predict  the  next  value  used  during  the  subsequent  perturbations.  An 
important  consideration  during  calculation  for  cases  near  the 
analytical  lift  limit  for  the  rotor  is  to  limit  the  extrapolated 
control  predictions  in  order  to  avoid  a  condition  too  far  into  stall 
(beyond  the  required  condition) .  The  two  Nawton-Raphson  methods 
which  produced  the  results  shown  here  do  not  include  specific  tests 
to  contend  with  predictions  which  overshoot  the  target  and  end  up 
too  far  into  stall.  This  would  be  a  problem  only  if  the  predicted 
controls  required  calculations  in  the  area  of  the  second  lift  rise 
at  very  high  angles  of  attack,  since  the  slopes  would  Indicate 
iterations  to  even  higher  angles  of  attack.  The  overshoot  at 
conditions  near  "maximum"  lift  is  most  critical  for  the  step-by-step 
conventional  method  since  convergence  is  accomplished  for  each 
individual  control  variable  (its  related  force,  direction,  or 
motion)  while  the  remaining  variables  are  held  constant.  This  means 
that,  if,  for  some  reason,  the  combination  of  controls  becomes 
unreasonable,  a  false  indication  of  lift  required  being  greater  than 
lift  available  will  result.  For  this  reason  some  checks  are 
required  which  result  in  a  restart  at  a  more  reasonable  value  for 
the  step-by-step  method. 

4.  Results 


4.1  Preliminary  Investigations 

Before  proceeding  to  the  results  for  a  production  airloads  program, 
it  is  interesting  to  compare  results  for  a  research  oriented 
response  problem  as  given  in  Reference  2.  Three  separate 
assumptions  can  be  established  for  a  comparison  of  periodic  shooting 
with  numerical  Integration  (solution  of  aquations  of  motion  until 
tranlents  decay).  Figure  8  Illustrates  the  boundaries  established 
when  these  assumptions  are  coupled  to  a  knowledge  of  the  stability 
of  problem.  For  the  sake  of  comparison,  it  is  assumed  that  each 
method  starts  with  a  first  goass  of  the  initial  conditions  (often 
zero),  having  an  error,  E  .  Each  method  moat  then  be  pursued  until 
a  desired  error,  E  is  reached.  If  is  also  assumed  that  the 
equations  are  nonlinear  so  chat  the  shooting  method  requires  several 
iterations.  For  the  controls  known  case.  Figure  8  shows  that  for 
10Z  damping,  such  as  is  typical  of  articulated  rotors  with  dampers, 
direct  numerical  integration  is  always  preferred,  even  for  only  one 
degree  of  freedom.  For  hlngeless  rotors,  however,  for  which  as 
little  as  1Z  damping  is  typical,  direct  integration  is  superior  only 
when  more  than  16  degrees  of  freedom  are  present.  For  damping  less 
chan  0.1Z,  as  is  typical  in  stability  work,  direct  integration  is 
generally  inferior  to  Che  present  method  of  periodic  shooting. 
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Figure  8.  Comparison  of  Numerical  Integration  and  Periodic  Shooting 


For  the  case  when  the  controls  are  not  known,  the  above 
comparison  must  be  modified  to  include  the  fact  that  controls  and 
initial  conditions  are  found  simultaneously  by  the  combined  Newron- 
Raphson  shooting  method,  but  are  found  sequentially  when  Newton- 
Raphson  is  coupled  with  direct  integration.  Figure  8  shows  for  the 
comparison  between  shooting  and  direct  integration  becomes  more 
favorable  for  periodic  shooting  when  Che  controls  must  be  found. 

For  typical  articulated  rotors  (damping  102)  shooting  is  superior  for 
less  than  10  degrees  of  freedom  and  for  typical  hingeless  rotors 
(damping  12)shooting  is  superior  for  less  than  100  degrees  of 
freedom.  Thus  there  is  a  great  potential  advantage  of  the  shooting 
method  over  numerical  integration  even  fcr  large  problems. 

Finally,  it  might  be  argued  that  the  potential  advantage  of 
direct  integration  would  increase  if  direct  integration  were  used 
with  only  an  estimated  set  of  partial  derivatives.  However,  as  seen 
in  Figure  8,  even  for  estimated  derivatives.  It  has  been  found  that 
there  is  still  a  favorable  trade-off  between  shooting  and 
integration. 

Vhus,  the  relative  advantage  of  shooting  is  enhanced  for 
systems  with  low  damping.  For  unstable  systems,  (damping  less  than 
0.0)  direct  integration  cannot  be  used  and  so,  periodic  shooting  (or 
some  other  method)  is  necessary. 
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In  order  co  provide  a  comparison  of  che  relative  efficiencies 
of  che  iteration  procedures  which  would  be  indicative  of  what  is 
required  co  undertake  a  complete  rotor  loads  and  performance 
analysis,  calculations  have  been  made  of  an  advanced  rotor  design 
for  which  wind  tunnel  data  is  available  (References  10,  11,  12). 

The  simulation  of  che  rectangular  planform  rotor  (baseline)  requires 
che  ability  co  include  three  airfoil  data  tables  with  interpolation 
between  adjacent  ones  co  represent  transition  section 
characteristics  (Figure  9) .  The  airfoil  data  cables  consist  of  data 
measured  in  a  transonic  wind  tunnel  at  Reynolds  numbers  which  are 
representative  of  che  full  scale  rotor  test. 
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Figure  9.  Blade  and  Tip  Geometry 


Figures  10,  11  and  12  illustrate  the  range  of  rotor  test 
conditions  and  che  number  of  points  which  ware  calculated  co  provide 
a  comparison.  The  actual  reported  test  point  values  (lift,  drag, 
side  force  and  shaft  angle)  are  used  as  the  calculated  trim  valves. 
This  not  only  exercises  the  trim  procedures  to  che  maximum  possible 
extent  for  this  analysis,  but  also  insures  che  calculation  of  the 
actual  rotor  condition  as  measured  In  che  test.  Although  the 
calculated  value  of  the  relative  rotor  power  does  not,  in  all  cases, 
correlate  well  with  the  test  value,  che  trends  are  quite 
representative  for  the  range  from  the  aucorotative  to  the  propulsive 
state  of  the  rotor. 


RELATIVE  ROIt>H  HOW  EH.  CH/S 


Correlation  can  be  obtained  through  the  adjustment  of  profile 
drag  and  inflow  velocity.  The  profile  drag  increment  would  account 
for  the  differences  between  the  wind  tunnel  airfoil  and  the  full 
scale  section  (small  imperfections).  The  inflow  velocity  varible 
adjustment  can  be  used  (References  13 ,  14)  to  adjust  the  3lope  of 
the  variation  of  relative  rotor  power  with  advance  ratio.  Through 
the  comparison  of  the  analysis  with  the  full  range  of  the  test 
results,  it  is  insured  that  the  trim  iteration  methods  are  exercised 
to  their  useful  limits.  The  test  results  represent  a  helicopter 
rotor  at  its  lift  and  propulsive  force  limits  (within  the 

power  required  limit  of  the  test  facility)  for  a  vide  range  of 
inflow  conditions. 

A  summary  of  the  rotor  revolutions  required  for  each  case  is 
shown  in  Figure  13  for  a  comparison  of  each  of  the  trim  iteration 
methods.  It  la  apparent  that  for  this  analysis,  the  Newt on- Rap h son 
and  Shooting  methods  are  superior  in  an  overall  reliability  and 
efficiency  sense.  However,  it  is  interesting  to  note  chat  the 
parallel  shooting  method  fails  to  converge  at  some  isolated  cases 
for  a  very  high  lift  condition,  where  Che  other  two  methods  are 
successful. 


Figure  13.  Summary  of  Convergence 


5.  Conclusions 


The  shooting  method  can  be  applied  to  an  existing  airloads 
program  with  a  moderate  amount  of  program  changes.  Many  airloads 
programs  are  set  up  to  remain  in  the  azimuth  loop  until  blade  motion 
transients  decay.  This  loop  must  be  interrupted  to  allow  only  one 
revolution  per  perturbation  and  the  resulting  partial  derivative 
then  added  to  the  matrix  for  the  Newton-Raphson  method. 

For  this  example,  parallel  shooting  is  superior  to  the 
conventional  method  for  about  30%  of  the  cases.  This  is  consistent 
with  earlier  estimates  for  a  system  with  1  Degree  of  Freedom  and 
14%  damping. 

The  use  of  an  approximate  Partial  Derivative  Matrix  is  not 
satisfactory  and  requires  an  average  of  2  -  8  times  as  many 
iterations  as  when  using  the  full  matrix. 

No  unusual  convergence  problems  were  encountered.  The  parallel 
shooting  method  and  the  conventional  method  generally  failed  to 
converge  for  the  same  cases  (about  6%  of  the  time) . 

The  sequential  application  of  shooting  is  generally  much  less 
efficient  than  the  parallel  application,  requiring  3  to  4  times  as 
many  rotor  revolutions. 

The  periodic  shooting  technique  can  be  successfully  applied  to 
a  rotor  airloads  program  which  includes  detailed  aerodynamics  and 
dynamic  stall,  to  calculate  the  full  range  of  performance  of  an 
advanced  technology  operational  helicopter  rotor. 
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